OSGeo / gdal

GDAL is an open source MIT licensed translator library for raster and vector geospatial data formats.
https://gdal.org
Other
4.81k stars 2.52k forks source link

MSGN driver is not able to open SEVIRI products in MSG1.5 Native format after a ROI has been applied #7915

Open mdbmdb74 opened 1 year ago

mdbmdb74 commented 1 year ago

Expected behaviour and actual behaviour.

GDAL 3.6.2. (also applies for 3.5.3, 3.7.0) is not able to open SEVIRI products in MSG1.5 Native format when a subset by a RegionOfInterest (ROI) has been applied: in fact, the EUMETSAT DataCenter allows users to require a ROI must be applied to the output products.

GDAL 3.3.2, instead, was able to correctly open those kind of files.

Steps to reproduce the problem.

The attached .zip file include the MSG1.5 Native file MSG3-SEVI-MSG15-0100-NA-20230607004509.766000000Z-20230607005759-4851554.nat as download from the EUMETSAT DataCentre. By means the DataCentre webapp a subset by RegionOfInterest has been applied to the output product: UL=29d32' 1.03"E - 48d13'54.00"N, LR=37d39'37.68"E - 43d33'23.21"N.

Using GDAL 3.5.3, 3.6.2 and 3.7.0 the following command

gdalinfo MSG3-SEVI-MSG15-0100-NA-20230607004509.766000000Z-20230607005759-4851554.nat

generates this error:

ERROR 1: Neither Whole Disk nor RSS - don't know how to handle.

Using GDAL 3.3.2 the same command generates the right output:

Driver: MSGN/EUMETSAT Archive native (.nat)
Files: /data/dt_data_for_test/HRSEVIRI/MSG3-SEVI-MSG15-0100-NA-20230607004509.766000000Z-20230607005759-4851554.nat
Size is 224, 113
Coordinate System is:
PROJCRS["Geostationary projection (MSG)",
    BASEGEOGCRS["MSG Ellipsoid",
        DATUM["MSG_DATUM",
            ELLIPSOID["MSG_SPHEROID",6356583.8,295.48806589701,
                LENGTHUNIT["metre",1,
                    ID["EPSG",9001]]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433,
                ID["EPSG",9122]]]],
    CONVERSION["unnamed",
        METHOD["Geostationary Satellite (Sweep Y)"],
        PARAMETER["Longitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Satellite Height",35785831,
            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]]]]
Data axis to CRS axis mapping: 1,2
Origin = (1944261.193359375000000,4362586.072753906250000)
Pixel Size = (3000.403076171875000,-3000.403076171875000)
Metadata:
  ch01_cal=-1.114656032994e+00 2.185600064695e-02
  ch02_cal=-1.465775717050e+00 2.874070033431e-02
  ch03_cal=-1.211260244250e+00 2.375020086765e-02
  ch04_cal=-1.865920103496e-01 3.658666869601e-03
  ch05_cal=-4.242236706827e-01 8.318111189856e-03
  ch06_cal=-1.969720376950e+00 3.862196817550e-02
  ch07_cal=-6.463960252982e+00 1.267443186859e-01
  ch08_cal=-5.302006445576e+00 1.039609106976e-01
  ch09_cal=-1.045681948659e+01 2.050356762077e-01
  ch10_cal=-1.133786848073e+01 2.223111466809e-01
  ch11_cal=-8.037951740768e+00 1.576068968778e-01
  Date/Time=20230607/00:57
  Origin=3310 1208
  Radiometric parameters format=offset slope
Corner Coordinates:
Upper Left  ( 1944261.193, 4362586.073) ( 29d32' 1.03"E, 48d13'54.00"N)
Lower Left  ( 1944261.193, 4023540.525) ( 26d11' 7.36"E, 42d43'49.71"N)
Upper Right ( 2616351.482, 4362586.073) ( 43d54'51.30"E, 49d30'31.53"N)
Lower Right ( 2616351.482, 4023540.525) ( 37d39'37.68"E, 43d33'23.21"N)
Center      ( 2280306.338, 4193063.299) ( 33d32'32.50"E, 45d48' 3.65"N)
Band 1 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 01
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 2 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 02
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 3 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 03
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 4 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 04
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 5 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 05
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 6 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 06
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 7 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 07
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 8 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 08
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 9 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 09
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 10 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 10
  Min=1.000 Max=1023.000 
  NoData Value=0
Band 11 Block=224x1 Type=UInt16, ColorInterp=Undefined
  Description = band 11
  Min=1.000 Max=1023.000 
  NoData Value=0

Operating system

Tested on both Rocky Linux 8.7 (Green Obsidian) and CentOS Linux 7 (Core)

GDAL version and provenance

GDAL 3.5.3, 3.6.2 from Anaconda and and 3.7.0 from conda-forge

Sample file as download form the EUMETSAT DataCentre

MSG3-SEVI-MSG15-0100-NA-20230607004509.766000000Z-20230607005759-4851554.nat.zip

Related issue

Maybe the new behaviour (error) has been introduced resolving the issue #1040.

rouault commented 1 year ago

@g8sqh The issue reported in this ticket seems to originate from your changes in #6244 / 3273d72713d . Will you be able to have a look ?

g8sqh commented 1 year ago

I will certainly look. Always wonder if anyone else uses these drivers - seems they do!

--David

On Fri, 09 Jun 2023 14:01:08 +0000 (UTC) Even Rouault @.***> wrote:

@g8sqh The issue reported in this ticket seems to originate from your changes in #6244 / 3273d72713d . Will you be able to have a look ?

-- Reply to this email directly or view it on GitHub: https://github.com/OSGeo/gdal/issues/7915#issuecomment-1584624609 You are receiving this because you were mentioned.

Message ID: @.***>

-- @.***>

mdbmdb74 commented 1 year ago

@g8sqh you might like to know that this driver is widely used in EUMETSAT's DataTailor webtool, a tool for which I am one of the developers.

rouault commented 1 year ago

you might like to know that this driver is widely used in EUMETSAT's DataTailor webtool, a tool for which I am one of the developers.

good to know. One issue with this driver is that we have absolutely zero regression tests for it. That would be a welcome addition if that driver must be relied upon. Given that test data files might be a bit big to be hosted in that repository, they should be hosted somewhere (potentially on http://download.osgeo.org/gdal/data/) and downloaded on the fly like done in other drivers (cf https://github.com/OSGeo/gdal/blob/master/autotest/gdrivers/nitf.py#L6226)

g8sqh commented 1 year ago

Oh ****

I receive EUMETCAST mainly for fun (a successor in a way to the VHF APT downlink) - and I hacked the driver to get the common cases working (0 deg, RSS, IODC and the hi-res pan channel on each) - and there were some small 'magic offsets' I needed to add (and never liked) to get good registration to ground truth. I just wanted to combine MSG with other sources in QGIS (and hence in GDAL)

Time to move up a level and make this a 'proper' driver if it's being used seriously. I will look at the regression test mechanism as suggested, but first I need to study how ROI is implemented. Certainly I only ever tested whole images, and particularly the split-window used on the pan channel may well implicitly assume size on image. (And even more need for regression - the window moves during the day sometimes...)

--David

On Fri, 09 Jun 2023 14:43:22 +0000 (UTC) mdbmdb74 @.***> wrote:

@g8sqh you might like to know that this driver is widely used in EUMETSAT's DataTailor webtool, a tool for which I am one of the developers.

-- Reply to this email directly or view it on GitHub: https://github.com/OSGeo/gdal/issues/7915#issuecomment-1584698650 You are receiving this because you were mentioned.

Message ID: @.***>

-- @.***>