Image co-registration with AROSICS

[Noza1212]

Dear Kvos.

I know Coastsat doesn´t correct the images and their correction is made by the data provider. But, I've noticed there´s an important gap or mismatch (up to 50 m in some images) between the images downloaded and an orthorectified image I got myself.

What would you suggest in that case? Is it ok if I rectify the images downloaded, delete the wrong ones and replace them?

Thank you in advance for your help.
Noza1212

[kvos]

hi @Noza1212 , this is a very good point. For Landsat images, the georeferencing accuracy is stored in the metadata (also in the output with the shorelines) and you can discard the images with insufficient georef accuracy. However, for Sentinel-2, which are more affected by these shifts, there isn't any georef error quantified.
You can georeference the tiffs before mapping the shorelines. This library is great for doing this https://pypi.org/project/arosics/, I have it on my todo list to implement this in coastsat one day.
let me know how you go as I'm really interested in this issue
cheers,
kilian

[Noza1212]

Thank you, I'll look into it.
Is there a manual or tutorial to learn how to use it?

[kvos]

yep, doco is here https://danschef.gitext-pages.gfz-potsdam.de/arosics/doc/ and there's also a paper on it https://www.mdpi.com/2072-4292/9/7/676

[Noza1212]

Hi @kvos

I've tried to use arosics in order to correct the images. It's been quite useful so far, I've got great results with Sentinel 2 images. Nonetheless, I couldn`t get a result with Landsat images. I think it is because of the way Coastsat downloads Landsat images.
Perhaps you can try it yourself when you have some free time.

[kvos]

sounds great, thanks for having a look! I can modify the way the images are downloaded if you let me know which bands/metadata is needed. Could you share a snippet of how you used the library to co-register the S2 images?

[Noza1212]

Well, as I undertand, basically, arosics has two classes: COREG and COREG_LOCAL. The difference between these two is the method applied.

COREG_LOCAL uses a grid (whose resolution needs to be defined first) in order to detect geometric shifts in a moving-window manner, so it works through the whole image. On the other hand, COREG focuses on a single coordinate position (which you can define), but its faster than COREG_LOCAL (of course, the speed of the process will depend on the resolution of the grid and the size of the image).

For example, this is an image that arosics creates for the user in order to visualize shifts. This was done using a Sentinel-2 image using COREG_LOCAL (you can see the grid as a dashed line in the image):

COREG basically needs the paths of the reference image and the target image as arguments, filling the params 'im_ref', 'im_tgt'. Additionally, you can define the path where you want to save the coregistered image filling the param 'path_out', or the format of the image with the param 'fmt_out' (GTIFF in case you want a GeoTiff as a result). You can also define the position coordinates of the window where you want the shifts detected with the param 'wp', and even the size of the window with the param 'ws'.

COREG_LOCAL needs the paths of the reference image, the target image and the resolution of the grid in order to work.

For example, you can use arosics like this:

from arosics import COREG
CR = COREG(im_ref, im_ tgt, wp= 523860, 9782890, ws= (128, 128), path_out= 'auto', fmt_out = 'GTIFF')

from arosics import COREG_LOCAL
CRL = COREG_LOCAL(im_ref, im_tgt, grid_res= 25, path_out= 'auto', fmt_out = 'GTIFF')

But, in order to get a result, afterwards, you need to calculate the spatial shifts. you can do it like this:
(If you used COREG) -> CR.calculate_spatial_shifts()
(If you used COREG_LOCAL)-> CRL.calculate_spatial_shifts()
Although it is the same function. it works kinda differently depending on the method.

In the end, you use the function correct_shifts to obtain your coregistered image in the format you defined.
(If you used COREG) -> CR.correct_shifts()
(If you used COREG_LOCAL)-> CRL.correct_shifts()

So finally, I've tried to made a synthesis of the arosics basics. As for your question about the bands needed, I am not an expert in that topic, so I don't know what to answer. However, I think the problem is there because when I read the Landsat images as array, the array is filled partially with '-inf', and when you try to calculate the spatial shifts, you get these error: RuntimeWarning: invalid value encountered in true_divide. That doesn´t happen with the Sentinel-2 images.

I hope this can be helpful for you. I'm really interested in solving this problem so if you have any other question, feel free to ask. I will keep trying to get a result.

[Noza1212]

Hi Kvos. I figured out the problem with Landsat images. It turns out that Landsat images come with missing data (-inf values) while Sentinel images do not.
So, you can solve the problem with the following line at the start of the code -> "np.seterr(divide='ignore', invalid='ignore')"

And with that, it should work just fine. One thing I noticed, though, was that with Landsat images the result of the co-registration process was not as good as the result I obtained with Sentinel 2 images. I think it's because Sentinel 2 images have a better resolution.

[kvos]

hi @Noza1212 , thanks very much for reporting this, it is super helpful! I only got to read everything in detail now.
Just some brainstorming here to see how to proceed with this enhancement, which I think will be a significant improvement!

• Does CORREG_LOCAL apply a different shift to each point in the grid? like correction for distortions in the image?
• Did you get better results with CORREG_LOCAL or CORREG?
• if using COREG we would need to find the pixels that are water automatically, which could be done using the coastsat classifier which has water and land classes.
• Finally, when adding this feature to the current toolbox, it could be done in pre-processing (co-registering the .tif before mapping the shorelines) or also in post-processing as the XY shifts could be applied to the shoreline coordinates directly without having to save an extra co-registered geotiff. I would be more keen for the post-processing option so that we can use the classifier to set the grid.

Please anyone who is interested in this feature feel free to participate to the discussion

[Noza1212]

Hi Kvos. Sorry for the delay in answering to you. I've been kinda busy these days.
In my short experience using Arosics what I can tell you is:

1. Yes. In fact, it calculates the shifts for every point in the grid and then, corrects the distortions according to the values of the shifts. Remember that, in order to get a good result with Arosics, the input images must already been thoroughly registered.

2. Well, for Coastsat, since you want to get the complete shoreline of any coast, it is not really useful if, for example, you get a good registration of the images at the north of the zone that you are studying and a deficient one at the south. COREG only focuses on a part of your image, you choose wich one (it could be the upper side, center, left center, lower center or wichever you want).
   However, COREG_LOCAL works on all the image, so I suggest you to use that function to enhance Coastsat.

3. I think I didn't understand what you meant with using COREG and finding water pixels. But I can tell you that the algorithm works o the land side of the image, because on water, it can't find any points to tie. I suggest, before downloading the images, try making the polygon in a way that it shows more land and less of the sea side in the image.

4. I think it's a great idea to apply it in the post-processing, since, as you say, you wouldn't need to save another image. In fact, Arosics gives you the option to get a table with the shifts values applied to each coordinate.

Finaly, I'm sorry if my english wasn't good enough. It's not my native language, so if I didn't explain myself good enough, please ask me again. I'm happy to help you with the enhacement of Coastsat, since it has been very useful to me.

[kvos]

thanks for the reply, was very clear. Sounds good for this enhancement, hopefully a first version of this enhancement can developed this year on the development branch, and we can start testing it. I'm in my last year of PhD so a bit tight on time but happy to contribute if anyone gets this started :)

[JanTiede]

Hi guys,
I also played around with arosics and worked with the corrected images. I don't have any data to validate the results. Kilian, is it possible for you to share the data of Narrabeen beach you have?

My understanding is that we would need to mask the water pixels and shoreline + ~50m perimeter to avoid the image being corrected for shifts in the shoreline, is that correct?

Cheers
Jan