I have written a program to perform subpixel image registration using mathematica and doing some of the more computationally complex parts of the code in CUDA for a nice speedboost. I have this all working fine and i can compute the misalignment between images to arbitrary fractions of pixels.
The next step for me is to take the misaligned picture and generate a new image from it which is perfectly aligned with my reference image. To do this i need to generate new values for each pixel using some sort of interpolation as they are effectively inbetween the pixels of the actual image. Sorry is this is worded badly…
It is my understanding that I should be able to do this pretty efficiently using the texture memory to hold the image as it can perform bilinear interpolation when i try to read values from floating point coordinates that are inbetween the actual pixels??
Assuming this is correct im just not sure how to do it. I get the feeling it should be pretty easy. As I am working from mathematica I can only load things into the global memory and im not sure where to proceed from here, do I make a kernel where I pass in the array of floating point values, then copy it to texture, then read back from the texture into a new array of floating point values? I am also unsure as to what will happen at the boundary conditions. I could probably just scrap the outermost pixel on every edge easily enough because I am working with pretty large images here.
Bilinear interpolation is indeed provided by CUDA’s texture support, however limited in accuracy (I think the interpolator is 8 bit quantized)
Consider using something more involved than bilinear interpolation.
The mathematically perfect way to do it would be a sinc filter, but that requires an infinite number of samples (in theory). There are other ways like bicubic interpolation which are still manageable computationally and produce results that are superior to bilinear interpolation. You might even find example implementations in CUDA on these forums…
Thanks for the info, the Lancsoz interpolation seems like it shouldn’t be too hard to implement although I’m unsure as to how badly any of these methods will affect the data. The images I’m registering are from a TEM and some are 32bit floating point if it will make a considerable difference.
I was originally just going to average the results from the 4 pixels that would effectively make up the new pixel weighted by the area within the new pixel if that makes any sense. I hadn’t really considered aliasing errors. My best bet I think is to try as many methods as I can and to test the results empirically.
Ideally the processing should be as fast as possible aswell as the rest of my program runs in a very short amount of time with image sizes upto around 4000*3000.