GFZ / arosics

AROSICS - Automated and Robust Open-Source Image Co-Registration Software
https://git.gfz-potsdam.de/danschef/arosics
Apache License 2.0
146 stars 26 forks source link

SSIM fails with window shape mismatch error #4

Closed adamjson64 closed 3 years ago

adamjson64 commented 3 years ago

Description

Attempting to coregister two rasters results in window size mismatch errors during SSIM calculation (off by one rounding errors) that result in the job failing. I was able to fix the problem on the dataset linked below by modifying the following code in the Deshifter class:

Before

xmin = find_nearest(self.out_grid[0], xmin, roundAlg='on', extrapolate=True, tolerance=x_tol)
ymin = find_nearest(self.out_grid[1], ymin, roundAlg='on', extrapolate=True, tolerance=y_tol)
xmax = find_nearest(self.out_grid[0], xmax, roundAlg='off', extrapolate=True, tolerance=x_tol)
ymax = find_nearest(self.out_grid[1], ymax, roundAlg='off', extrapolate=True, tolerance=y_tol)

After

xmin = find_nearest(self.out_grid[0], xmin, roundAlg='off', extrapolate=True, tolerance=x_tol)
ymin = find_nearest(self.out_grid[1], ymin, roundAlg='off', extrapolate=True, tolerance=y_tol)
xmax = find_nearest(self.out_grid[0], xmax, roundAlg='off', extrapolate=True, tolerance=x_tol)
ymax = find_nearest(self.out_grid[1], ymax, roundAlg='off', extrapolate=True, tolerance=y_tol)

Dataset https://drive.google.com/file/d/1Ftc9tLTAz3H1dHVx-GHr0KgUTaSzW_cw/view?usp=sharing

What I Did

im_reference = 'rgb.tif'
im_target = 'red.tif'
kwargs = {
    'grid_res': 512,
    'window_size': (2056, 2056),
    'path_out': '/tmp/coregistration/output.tif',
    'projectDir': None,
    'q': False,
    'v': True,
    'fmt_out': 'GTIFF',
    # 'tieP_filter_level': 3,
    # 'max_points': 100000,
    'min_reliability': 60,
    'max_shift': 128,
    'max_iter': 100,
    'match_gsd': False,
    'out_crea_options': ['PHOTOMETRIC=MINISBLACK',
                         'COMPRESS=LZW',
                         'ALPHA=NO']
}
kwargs['nodata'] = (-10000, -10000)
CRL = COREG_LOCAL(im_reference, im_target, **kwargs)
CRL.correct_shifts()
  File "process.py", line 125, in process
    CRL.correct_shifts()
  File "/usr/local/lib/python3.6/dist-packages/arosics/CoReg_local.py", line 754, in correct_shifts
    self.calculate_spatial_shifts()
  File "/usr/local/lib/python3.6/dist-packages/arosics/CoReg_local.py", line 437, in calculate_spatial_shifts
    self._tiepoint_grid.get_CoRegPoints_table()
  File "/usr/local/lib/python3.6/dist-packages/arosics/Tie_Point_Grid.py", line 400, in get_CoRegPoints_table
    results = results.get()
  File "/usr/lib/python3.6/multiprocessing/pool.py", line 644, in get
    raise self._value
ValueError: Input images must have the same dimensions.
danschef commented 3 years ago

Thanks for filing the issue. I will take a look at it as soon as possible.

danschef commented 3 years ago

Is there any reason why you use a window_size of 2056x2056?

The default is at 256x256 and I would not recommend to use window sizes larger than 512x512. The coregistration accuracy won´t increase but processing takes much longer.

Apart from that, AROSICS will automatically reduce the window size to avoid including no-data areas or areas outside of the image overlap. That means, that the final window size used for all the computed tie points varies a lot.

adamjson64 commented 3 years ago

The thought was that since its high res drone imagery the wider search area would ensure there are enough features to allow for an accurate match and empirical evidence seemed to back that theory up.

Is my logic flawed?

danschef commented 3 years ago

Ok thanks, that may make sense. I did not have a look at the data, yet.

adamjson64 commented 3 years ago

I'm not actually sure that my "fix" fixes this in all cases. I'm digging into it now..

adamjson64 commented 3 years ago

Is it possible using gdalwarp to migrate from one CRS to another could play a role here? Maybe thats a long shot but im grasping at this point. The only difference between the dataset with the issue and all the other datasets I've successfully processed is that this dataset was warped to a different CRS

danschef commented 3 years ago

I increased the tolerance used for snapping to the coordinate grid for now in https://github.com/GFZ/arosics/commit/242bdb5fc46ea5c59bae1fce0c78493b14d0e2c3. That fixes the issue but I am not completely sure if that is the optimal solution. If you have a better solution, let me know. I think it is related to inaccuracies when transforming float coordinates from one CRS to another. This may cause errors in the range of a small fraction of a pixel which then lets the coordinate snap to a different coordinate on the grid. Using a tolerance one 2000th of a pixel instead of one 10000th fixed it.

You may update to AROSICS version 1.4.2 soon. But make sure to update it with conda or by running pip install (only works if the requirements are already there) because there are also quite some changes in my upstream packages.

Thanks for providing the test data. I could easily reproduce the issue.

Be aware that you passed a wrong no-data value for the reference image. It is 0 and not -10000. This caused AROSICS to compute a wrong nodata mask which may result in some inaccurate tie points at the edge of the overlap area. I added a validation check that now rejects out-of-range nodata values (data type of your reference image us uint8, so -10000 won´t work).

Also, it makes not much sense to increase the number iterations for finding a match to 100. In moss cases, there is only one iteration needed - using 100 simply slows down processing without much of improvement.

To really improve the co-registration result from now, I think you have to go back one step to the mosaicking of your drone flight paths. The reference image has quite some spatial artifacts, especially in the lower left corner - there are areas of the ground surface that are missing in the image, I guess due to errors in mosaicking. These artifacts make it hard to find useful matches there.

Here is a screenshot of the current output: image

Also take a look to the 3-step filtering techniques. I may make sense to adjust some thresholds there to improve the results.

danschef commented 3 years ago

Closing this for now, please re-open if there are any issues.

adamjson64 commented 3 years ago

Thanks so much for fixing this!

I actually left the incorrect nodata value as a hack because only 57% of the reference image was said to be overlapping and I was afraid that would cause more harm than the regions outside the raster area having tie points. I think this is probably due to a few small no data areas inside the raster. Do you have any suggestions on how to deal with this? Should I attempt to fill those nodata areas?

danschef commented 3 years ago

Check out how this works with the new version when you provide nodata=(0, -10000). However, the new releases are still in progress but they should be available in the course of the day.

adamjson64 commented 3 years ago

Yeah I made nodata the change on my end - thanks. Just a heads up - auto detection is broken for the dataset associated with this ticket. It finds 0.0 as the nodata value but throws an exception presumably because its a float instead of an int.

When do you expect version 1.4.2 to be available in conda-forge?

danschef commented 3 years ago

This is already fixed in 1.4.3. Should be available in conda-forge now.