Infoseite zur "schnellen" normierten Kreuzkorrelation im Rahmen der Vorlesung BV1 im WiSe 2008/09

Materialien:

• Der Vortrag vom 25.11.2008 zu dem Thema schnelle normierte Kreuzkorrelation.
• Das Paper zur schnellen Berechnung der normierten Kreuzkorrelation, auf dem der Vortrag aufbaut: J. P. Lewis - Fast Normalized Cross Correlation, 1995

...speziell für die schnelle Bildverarbeitung unter Mac OS X:

• Einen Überblick über das Accelerate-Framworks kann man bei Apple unter Taking Advantage of the Accelerate Framework finden.
• Ein gutes Tutorial zu dem Accelerate Framework und den Vorgehensweisen (in mehreren Teilen) findet man bei MacResearch:
  • Teil 1: Einführung
  • Teil 2: Benutzung der Funktionen (1)
  • Teil 3: Benutzung der Funktionen (2)
• Die komplette Dokumentation zu der 2D-Fast Fourier Transform in dem Accelerate-Framework stellt Apple unter der vDSP Two-Dimensional Fast Fourier Transforms Reference zur Verfügung.

Benutzung der FFT des Accelerate Frameworks unter Mac OS X (Objective-C)

Um ein Bild an die FFT zu übergeben, muss es zunächst folgendes erfüllen:

• Es sollte in der richtigen Größe vorliegen: Im Idealfall: Länge = Breite = 2^n
  Sollte dies nicht der Fall sein, muss das Bild evt. vorher in ein größeres kopiert werden.
• Es muss ein einzelnes C-Array vom Typ DSPComplex sein.
  Dabei sollte es zeilenweise hintereinander (siehe Grafik) abgelegt sein.

0	1	2	3	4
1	0	1	2	3
2	4	5	6	7
3	8	9	10	11
4	12	13	14	15

Abbildung 1: Ein Bild mit Koordinatenachse (fett) und Scanadresse der einzelnen Pixel (Offset)

Die Umwandlung eines eingelesenen Bildes (genauer: eines CGImages) kann man hier nachlesen:
Technical Q&A QA1509 - Getting the pixel data from a CGImage object

Liegt das Bild als 1-Kanal-Repräsentation im geforderten Format vor, so kann man es an die FFT übergeben

#import <Cocoa/Cocoa.h>
#import <QuartzCore/QuartzCore.h>
#import <Accelerate/Accelerate.h>

int fft_forward(DSPComplex * img, int width, int height){
	return fft(img, width, height, FFT_FORWARD);
}

int fft_backward(DSPComplex * img, int width, int height){
	return fft(img, width, height, FFT_INVERSE);
}

int fft(DSPComplex * data, int width, int height, FFTDirection dir){
	
	/* Total area or volume covered by the dimensions */
	int dims = width*height;
	
	/* Get log (base 2) of each dimension */
	int logX = log2((double)width);
	int logY = log2((double)height);

	
	/* Setup and allocate some memory for performing the FFTs */
	DSPSplitComplex input;
	
	input.realp = (float*)malloc(dims*sizeof(float));
	input.imagp = (float*)malloc(dims*sizeof(float));
	
	/* Setup the weights (twiddle factors) for performing the FFT */
	FFTSetup fft_weights;
	fft_weights = create_fftsetup(MIN(logX,logY),kFFTRadix2);
	
	/* Split the complex (interleaved) data into two arrays */
	ctoz(data,2,&input,1,dims);
	
	/* Perform the FFT */
	fft2d_zip(fft_weights,&input,1,0,logX,logY,dir);
	/* Merge the separate data (real and imaginary) parts into a single array */
	ztoc(&input,1,data,2,dims);
	
	//Normalize energy 
	if(dir == FFT_INVERSE){
		int i;
		for(i=0;i<dims;i++){
			data[i].real=data[i].real/dims;
			data[i].imag=data[i].imag/dims;
		}
	}
	/* Free up the allocations */
	free(input.realp);
	free(input.imagp);
	
	destroy_fftsetup(fft_weights);
	
	return 0;
}

Download:

Selbstverständlich kann der Quelltext (deutlich umfangreicher als oben) hier heruntergeladen werden. Es wird ein Bild eingelesen, automatisch vergrößert und im Anschluß wieder abgespeichert. Dabei habe ich für die Bildrepräsentation allerdings Objekte vom Typ NSImage gewählt.