Co-registration of PlanetScope Imagery

I am using the AROSICS Python package to co-register PlanetScope images, see https://github.com/GFZ/arosics and the associated journal paper. I don't know whether this can be done with GEE.

@b.g.ruessink very useful, thank you for sharing!

In python, I have used the reproject_match command in rioxarray: rioxarray package — rioxarray 0.13.3 documentation (corteva.github.io). It appears to be based on the rasterio warp module (rasterio.warp module), which appears to be essentially the GDAL Warp tool (from the github page I think it just uses the GDAL Warp tool through python). Also not sure how it interfaces with GEE.

@tmyeoman Hi,

Depending on what you are comfortable with you can also using the orders API raster tool kit which also include a co-registration tool which you can simply add to your JSON payload. 

https://developers.planet.com/apis/orders/tools/#coregister

Thanks @Pooja Pandey .  Do you have any documentation or examples of how this tool might be applied to larger areas of interest?  The reason I didn’t use the co-registration tool in my order was because the page you linked suggested it was for small areas only.

From the page -  “This tool is designed to support coregistration of small areas of interest – contained within a single scene – and works best with high geographic overlap between scenes in the times series.”

@tmyeoman how large of an area are we talking as one superdove scene is around 600 sqkm? Thanks! If it is bigger then that, then i would revert to using external methods. 

@b.g.ruessinkhave you tried this with 1B_MS imagery, and do you need to do any grayscale reduction first?

We’re trying to do some stereo matching with pairs of PSScenes and the initial registration based on RPC alone is quite bad. For example I’ve had trouble getting enough SIFT matches for alignment over a city area (which seems odd - visually there’s plenty of uniqueness to match), and then the matchers crash out because there aren’t enough valid tiepoints.

(arosics) josh@:~$ arosics global -max_iter 10 -max_shift 128  files/PSScene/20201124_031355_67_1064/basic_analytic_udm2/20201124_031355_67_1064_1B_AnalyticMS.tif files/PSScene/20201126_031738_79_105d/basic_analytic_udm2/20201126_031738_79_105d_1B_AnalyticMS.tif
Calculating footprint polygon and actual data corner coordinates for reference image...
Bounding box of calculated footprint for reference image:
	(0.0, -4400.0, 6600.0, 0.0)
Calculating footprint polygon and actual data corner coordinates for image to be shifted...
Bounding box of calculated footprint for image to be shifted:
	(0.0, -4400.0, 6600.0, 0.0)
Matching window position (X,Y): 3300.0/-2200.0
No clear match found yet. Jumping to iteration 2...
input shifts:  -10 0
No clear match found yet. Jumping to iteration 3...
input shifts:  -30 0
No clear match found yet. Jumping to iteration 4...
input shifts:  -8 0
No clear match found yet. Jumping to iteration 5...
input shifts:  -32 0
No clear match found yet. Jumping to iteration 6...
input shifts:  -6 0
No clear match found yet. Jumping to iteration 7...
input shifts:  -34 0
No clear match found yet. Jumping to iteration 8...
input shifts:  -9 0
No clear match found yet. Jumping to iteration 9...
input shifts:  -31 0
No clear match found yet. Jumping to iteration 10...
input shifts:  -7 0
Traceback (most recent call last):
  File "/opt/conda/envs/arosics/bin/arosics", line 11, in <module>
    sys.exit(main())
             ^^^^^^
  File "/opt/conda/envs/arosics/lib/python3.11/site-packages/arosics/arosics_cli.py", line 363, in main
    parsed_args.func(parsed_args)
  File "/opt/conda/envs/arosics/lib/python3.11/site-packages/arosics/arosics_cli.py", line 69, in run_global_coreg
    COREG_obj.correct_shifts()
  File "/opt/conda/envs/arosics/lib/python3.11/site-packages/arosics/CoReg.py", line 1721, in correct_shifts
    DS = DESHIFTER(self.shift, self.coreg_info,
                               ^^^^^^^^^^^^^^^
  File "/opt/conda/envs/arosics/lib/python3.11/site-packages/arosics/CoReg.py", line 1661, in coreg_info
    self.calculate_spatial_shifts()
  File "/opt/conda/envs/arosics/lib/python3.11/site-packages/arosics/CoReg.py", line 1572, in calculate_spatial_shifts
    self._handle_error(RuntimeError('No match found in the given window.'))
  File "/opt/conda/envs/arosics/lib/python3.11/site-packages/arosics/CoReg.py", line 473, in _handle_error
    raise error
RuntimeError: No match found in the given window.

@joshvm I am working on 3B_AnalyticMS_SR imagery, for a relatively small study area of coastal dunes (a few sq km). In Arosics I am using the global co-registration (so not local) on the red band, and for most of images in my data set there are no issues. Occasionally, Arosics fails but I am a bit in the dark on why that is. I found that some of these images can be co-registered if I adjust the window size, for example, ws=(128,128) rather than ws=(256,256). Have you checked how the size of your search window, relative to building size, influences the performance of Arosics? I noted that you use a max_shift of 128. That is really large compared to the default of 5. Why such a large value? And, are there any issues with building-induced shadows that may confuse Arosics? My imagery sometimes has quite some topography-induced shadows (rugged terrain), but Arosics doesn't seem to be troubled by that.

Yeah this is a full tile pair over a coastal city. The default settings failed immediately so I wasn't sure if the shift was too small (hence setting to an enormous value). Though it seems you also have to set max iter as well, as they're dependent - eg 5 iters will only search a very small range regardless (as in the log, it doesn't search that widely).

There are some shadows but i really doubt it. qgis managed (and there are no options at all) so I don't think it should need much tweaking.

I'll try global/local, though global failed as above.Thanks!