# Ex14

Author: Jan Beeck

# Phase correlation with polar coordinates

Definition: Using phase correlation the angle can be easily found out.

Steps:
• Both images must be the same size.
• FFT from the images.
• Polar transform.
• FFT from the polar spectrum.
• Phase coorelation from the polar spectrum. (angle)
• Inverse FFT
1 def ictransformacaopolar(f,m,n):
2
3     m,n = f.shape
4     Rx = floor(m/2)
5     Ry = floor(n/2)
6
7     f[Rx,:] = (f[Rx-1,:]+f[Rx+1,:])/2
8     f[:,Ry] = (f[:,Ry-1]+f[:,Ry+1])/2
9
10     b = min(Rx,Ry)/m
11     a = pi/n
12
13     xl,yl = indices((m,n))
14
15     XI = Rx + (b*yl)*cos(a*xl)
16     YI = Ry + (b*yl)*sin(a*xl)
17
18     g = iainterpolclosest(f, (XI,YI))
19
20     return g
1 from ia636 import *
2 from numpy import *
3
4 def phase_correlation_p(f,g): #both images have to be the same
5     F = fft.fftn(f.astype(float))
6     G = fft.fftn(f.astype(float))
7
8     pF = ictransformacaopolar(F,F.shape[0],F.shape[1])
9     pG = ictransformacaopolar(G,G.shape[0],G.shape[1])
10
11     pFF = fft.fftn(pF)
12     pGG = fft.fftn(pG)
13
14     crossproduct = (pFF*pGG.conj())/(abs(pFF*pGG.conj()))
15     peak = fft.ifftn(crossproduct).real
16
17     #p, pos = max(ravel(peak)), argmax(ravel(peak))
18     #(row, col) = iaind2sub(peak.shape, pos)
19
20     return peak #=angle
1 from numpy import *
2
5 T = array([[1,0,0],[0,1,0],[0,0,1]], 'd')
7 Rotation = array([[cos(angle),-sin(angle),0],[sin(angle),cos(angle),0],[0,0,1]],'d')
8 T = dot(T,Rotation)
9 g = iaffine(f,T)
11
12 #Polar coordinates
13 angle = phase_correlation_p(f,g)
17 #print "peaks:", angle
18 maxvalue = unravel_index(argmax(angle), angle.shape)
19 ans = mod((float(maxvalue[1])/angle.shape[1])*2*pi,pi)
20 print "angle:", ans
angle: 0.0

# Volume Registration using the 3-D pseudo-polar Fourier

## Abstract

This paper introduces an algorithm for the registration of rotated and translated volumes using the 3-D pseudo-polar Fourier transform, which accurately computes the Fourier transform of the registered volumes on a near-spherical 3-D domain without using interpolation. We propose a three-step procedure. The first estimates the rotation axis. The second computes the planar rotation relative to the rotation axis, and the third recovers the translational displacement. The rotation estimation is based on Euler’s theorem, which allows to represent a 3-D rotation as a planar rotation around a 3-D rotation axis. This axis is accurately recovered by the 3-D pseudo-polar Fourier transform using radial integrations. The residual planar rotation is computed by an extension of the angular difference function to cylindrical motion. Experimental results show that the algorithm is accurate and robust to noise.

## Introduction

• There are several approaches towards 3-D registration: feature or intensity.

• Intensity based schemes align the volumes by using the intensities of the voxels. It also include optimizitation and Fourier based schemes.
• Frecuency domain: These methods use the properties of the Fourier transform to separetely estimate rotation and translation. This reduces a problem with six degrees of freedom into two problems with three degrees of freedom.

• First step: recovers the rotation axis.
• Second step: recovers the rotation angle.
• Third step: recovers the translation parameters using phase correlation.
• The 3-D spherical Fourier transform: the density volumes are aligned by computing the magnitude of their polar Fourier transform. As the volumes are given on Cartesian grids, the polar Fourier transform is interpolated from the Cartesian 3D Fourier transform coefficients. Thus, rotations are reduced to translations in the spherical coordinate system that are recovered by applying phase correlation. The residual translation is also estimated by phase correlation.

• A Fourier based approach: which does not require interpolation in the frequency domain. It is based on the 2-D pseudo-polar Fourier transform (PPFT2D) and 3-D pseudo-polar Fourier transform (PPFT3D).

• Compute the Discrete Fourier Transform (DFT) on non-Cartesian grids. Allowing a fast and algebraically accurate registration, which draws on Euler’s theorem to estimate the 3-D rotation.

• The algorithm has three steps:
• First: the 3-D pseudo-polar Fourier transform is used to recover the rotation axis.
• Second: the rotation around the axis is estimated using a pseudo-cylindrical representation computed with the 2-D pseudo-polar Fourier transform.
• Finally: the translation is computed by using phase correlation.

Note: The 3-D pseudo-polar Fourier transform (PPFT3D) evaluates the 3-D Fourier transform of a volume on the 3-D pseudo-polar grid (called P), which approximates the 3-D spherical grid.