kvos / CoastSat

Global shoreline mapping tool from satellite imagery
http://coastsat.space
GNU General Public License v3.0
679 stars 249 forks source link

Image co-registration with AROSICS #148

Open Noza1212 opened 3 years ago

Noza1212 commented 3 years ago

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 commented 3 years ago

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 commented 3 years ago

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

kvos commented 3 years ago

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 commented 3 years ago

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 commented 3 years ago

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 commented 3 years ago

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): Figure_1

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(imref, 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()

Figure_2

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 commented 3 years ago

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 commented 3 years ago

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!

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

Noza1212 commented 3 years ago

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 commented 3 years ago

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 commented 3 years ago

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

shannonmgbrown commented 2 years ago

Hi everyone, I am also trying to incorporate Arosics into the CoastSat workflow. Real quick first, I just want to say how awesome this tool is Killian. In trying to use arosics, I'm running into a couple hurdles that it sounds like you all might be experiencing too.

@Noza1212 , did you have any problem with CoastSat drawing a shoreline on the coregistered S2 image? I have successfully coregistered the 10m S2 tifs to a L8 image, but when I do the extraction step on CoastSat there is no shoreline drawn. The classification and MNDVI panels look good though. Did you experience this?

kvos commented 2 years ago

hi @JanTiede and @shannonmgbrown , cool that you have been using AROSICS, that's really exciting. Sorry for the late reply from my end. Yes I will make a new branch with some in situ data to validate your results (maybe Narrabeen to start with). Would be really keen to hear about your experiences with AROSICS and let's try and validate some of the results. Will send the new branch soon. Cheers, Kilian

kvos commented 2 years ago

@JanTiede , @shannonmgbrown , I created a new validation branch which has the in situ data for Narrabeen and a script validation_Narrabeen.py which performs the shoreline mapping (with automatic QA) + comparison with in situ surveys. check it out, it can be useful for validating the improvements with AROSICS coreg

kvos commented 2 years ago

please check out the new validation branch, I just updated it with a full jupyter notebook example of comparison against in situ data at Narrabeen https://github.com/kvos/CoastSat/blob/validation/validation_Narrabeen.ipynb

@JanTiede @shannonmgbrown

gianfrancodp commented 6 months ago

I just want to give a small contribution to think on this issue. I make a little test wit Arosics v.1.10.2. with two S2 scenes delayed of a half hour. In case of a time-series scenes I think that there could be linkage errors between the images which are perpetuated or even cancel out the final goals of analyzing the variation of the coastline.

I found another approach useful for CoastSat.

Andres Perez developed an ML tool for coregistration of a small time-series scenes: https://link.springer.com/article/10.1007/s00521-023-09242-0

repository is here: https://github.com/afperezm/multi-temporal-coregistration