Add_coastlines.sh removes most of geotiff headers in tif destination file

[ve2fxl]

Hi, I'm using add_coastlines.sh like this in my scripts:
add_coastlines.sh --add-coastlines --coastlines-resolution=i --coastlines-width 0.5 --add-borders --borders-resolution=i --borders-outline='red' --borders-width 0.5 *.tif -o *.tif -v

I even tested to output the destination file in another folder but the result is the same.
Geotiff infos before applying the coastlines:
Driver: GTiff/GeoTIFF
Files: GOES-16_ABI_RadF_C02_20230720_180021_GOES-East.tif
Size is 21696, 21696
Coordinate System is:
PROJCRS["unknown",
BASEGEOGCRS["unknown",
DATUM["unknown",
ELLIPSOID["GRS 1980",6378137,298.257222096042,
LENGTHUNIT["metre",1,
ID["EPSG",9001]]]],
PRIMEM["Greenwich",0,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8901]]],
CONVERSION["unknown",
METHOD["Geostationary Satellite (Sweep X)"],
PARAMETER["Longitude of natural origin",-75,
ANGLEUNIT["degree",0.0174532925199433],
ID["EPSG",8802]],
PARAMETER["Satellite Height",35786023,
LENGTHUNIT["metre",1,
ID["EPSG",9001]]],
PARAMETER["False easting",0,
LENGTHUNIT["metre",1],
ID["EPSG",8806]],
PARAMETER["False northing",0,
LENGTHUNIT["metre",1],
ID["EPSG",8807]]],
CS[Cartesian,2],
AXIS["easting",east,
ORDER[1],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]],
AXIS["northing",north,
ORDER[2],
LENGTHUNIT["metre",1,
ID["EPSG",9001]]],
REMARK["PROJ CRS string: +proj=geos +sweep=x +lon_0=-75 +h=35786023 +x_0=0 +y_0=0 +ellps=GRS80 +units=m +no_defs"]]
Data axis to CRS axis mapping: 1,2
Origin = (-5434894.885056000202894,5434894.885056000202894)
Pixel Size = (501.004322000000002,-501.004322000000002)
Metadata:
AREA_OR_POINT=Area
offset=nan
scale=nan
TIFFTAG_DATETIME=2023:07:20 18:00:21
Image Structure Metadata:
COMPRESSION=LZW
INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left (-5434894.885, 5434894.885) ERROR 1: Point outside of projection domain

Lower Left (-5434894.885,-5434894.885) ERROR 1: Point outside of projection domain

Upper Right ( 5434894.885, 5434894.885) ERROR 1: Point outside of projection domain

Lower Right ( 5434894.885,-5434894.885) ERROR 1: Point outside of projection domain

Center ( 0.0000000, 0.0000000) ( 75d 0' 0.00"W, 0d 0' 0.01"N)
Band 1 Block=256x256 Type=Byte, ColorInterp=Gray
Mask Flags: PER_DATASET ALPHA
Band 2 Block=256x256 Type=Byte, ColorInterp=Alpha

AND after applying the coastlines:
Driver: GTiff/GeoTIFF
Files: test.tif
Size is 21696, 21696
Image Structure Metadata:
COMPRESSION=LZW
INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left ( 0.0, 0.0)
Lower Left ( 0.0,21696.0)
Upper Right (21696.0, 0.0)
Lower Right (21696.0,21696.0)
Center (10848.0,10848.0)
Band 1 Block=21696x1 Type=Byte, ColorInterp=Red
Mask Flags: PER_DATASET ALPHA
Band 2 Block=21696x1 Type=Byte, ColorInterp=Green
Mask Flags: PER_DATASET ALPHA
Band 3 Block=21696x1 Type=Byte, ColorInterp=Blue
Mask Flags: PER_DATASET ALPHA
Band 4 Block=21696x1 Type=Byte, ColorInterp=Alpha

Is it a known issue and it's by design or a bug?

Is there a solution?

Thanks,
Luc Fontaine
VE2FXL

[davidh]

Hi,

The "add_coastlines.sh" script can only output PNG images. Additionally, the last released version may have had a bug where it didn't properly compare the number of output filenames to the number of input filenames and it would result in an error when you did it correctly. Regardless, I think the output files you are making are actually PNG images and the gdalinfo output is actually finding any geolocation information. You should use a ".png" extension for your output filenames and should always save the output to a new filename and not overwrite the input files. By default the script will take the input tif name and change the extension to `.png`. If you have any error messages from running "add_coastlines.sh", please let me know.

I hope this helps, but let me know if I misunderstood something.

Dave

[ve2fxl]

Hi Dave,

I open the resulting file with MediaInfo and it told that it was a tif file in fact. When I renamed the file to .png, MediaInfo told again it was a tif file and that the file extension should be .tif or .tiff
Why I switched from PNG to TIF0? Because the add_coastlines.sh script is quite faster outputing TIF than PNG and I thought also that I would be able to play with geotiff after the coastlines were added but that's not the case
Why I overwrite the input file with the output file? Because it just works ok, surprising but yes it works and that doesn't take more space (or maybe it does but hidden to me).

So there's no way that add_coastlines.sh will keep the geotiff infos in a future version?

Thanks,
Luc

[davidh]

Sorry for the confusion. It wasn't clear from your original message what wasn't working. So you are saying your main issue is that the geolocation information is not passed to the output file?

This is expected. My guess for why the output file is still considered a TIFF is that the library being used to write PNG/JPEG images is also capable of writing regular (not geo) TIFF images. However, this is not something we ever intended to support so the obvious question of "where did the geolocation information go" was not something I had considered handling. I've always assumed or thought that users should not be writing coastlines to geotiffs as the tools that geotiffs are usually used with allow you to add coastlines (Web Mapping Services, QGIS, ArcGIS, etc).

With the current version of the code it is not possible to get the geotiff metadata passed through to the output. You are the first to ask about this functionality, could you provide details about your use case? Note that even though you are providing the same filename for output, it is not updating the input image file but overwriting it (deleting and creating a new file). There may be a command line way of copying the geo metadata from the input file to the (separate) output file using GDAL tools but I've never done it before. If you can provide more details on your use case, I'll see how much work it would be to add this functionality.

[ve2fxl]

I have 2 cases I wanna work on.
Overlays here are all applied from 1 geo2grid server so when I want to do the next 2 things afterward, if the geotiff headers aren't available anymore, that won't work.
1-Cropping overlapping sections of overlapping mesos with gdal
2-Tracking center/eye of Tropical Storms/Hurricanes from NHC Tropical Cyclones Advisories (or better sources if available)

If not possible, I'll have to add the coastlines only after doing 1 or 2 or for any animation/still I make.

Thanks
Luc

[davidh]

Interesting. In the past I know users have done similar things as part of their Geo2Grid (or Satpy) processing, but I think in those cases it was their single main goal and not one of many possible outputs they wanted to make. Right now your only option is to swap the order of your processing or find a way to copy the geotiff tags from the source geotiff to the output geotiff after "add_coastlines.sh" is finished (GDAL maybe?). I've made a GitHub issue to track the progress of this feature though as I've discovered a way to copy the geotiff tags in Python in a seemingly simple way:

https://github.com/ssec/polar2grid/issues/627

[ve2fxl]

Thanks David,

I'll take a look!

Luc