Differences in projecting satellite images with cartopy and pyresample

Solution 1:

1. pyresample is much faster than cartopy (depending on the size of the output image up to a factor of 10)

The cartopy reprojection functionality hasn't been optimized in any way, and although it is using the scipy ckdtree functionality under the hood, the algorithm itself is written in Python. I seem to remember that a quick win was to use https://pypi.python.org/pypi/kdtree which from memory gave quite a reasonable speedup with little work, cartopy.img_transform would be the place where changes would be needed.

Cartopy's re-projection functionality is probably also paying the cost of being very general - you can provide an image in any projection, and it will put it into any other projection, dealing with discontinuities and tears without a problem. It would be really cool to hook into pyresample's functionality though (and GDAL's for that matter) to give users the opportunity to speed up the reprojection in certain cases.

2. The output images differ: The results using pyresample show a stronger contrast.

Looks like you're creating a matplotlib figure to resample the image and using mpl's savefig functionality. It is possible that this process is causing the contrast to be lost. I'd advise just using cartopy's reprojection functionality without adding an image to a figure and saving the figure (example at the end).

3. If the multiprocessing library is used to do the projection in parallel, the cartopy version crashes with the following error message:

This really surprised me as there is no C code in cartopy which is doing the reprojecting. Therefore you've either found a bug with scipy, or more likely you are hitting a problem with numpy/matplotlib (google brings up a few results with your exception and matplotlib and/or numpy, e.g. https://github.com/numpy/numpy/issues/1270).

Ok, so here is how I would do the reprojection without using matplotlib at all:

import cartopy.crs as ccrs
from cartopy.img_transform import warp_array
import numpy as np    
import PIL.Image


# I've downloaded the file from https://github.com/jmozmoz/cloudmap/blob/78923d15ad906eaa6d1dcab168a6364643d3fc94/debug/2014_8_7_1800_GOES15_4_S1.jpeg
# and clipped the image. 
fname = '2014_8_7_1800_GOES15_4_S1.jpeg'
img = PIL.Image.open(fname)

result_array, extent = warp_array(np.array(img),
                                  source_proj=ccrs.Geostationary(),
                                  target_proj=ccrs.PlateCarree(),
                                  target_res=(4000, 2000))

result = PIL.Image.fromarray(result_array)
result.save('reprojected.jpeg')

With the resulting image (eventually) looking something like:

There are some real possibilities for some optimisations with this functionality - quite a large amount of work is done creating the kdtree in the first place (which could potentially be cached) and another large chunk of the work is computing the indices from the original image (again, caches very well) which would essentially reduce and repeat reprojections to an numpy indexing problem.

If you want to look into the performance possibilities or the contrast issue (which I'm uncertain whether my solution fixes or not) please feel free to open up an issue on the github repo and we can talk through some of the options.

Thanks for asking, and HTH!