<div dir="ltr">Thanks! Educational as ever, this closes out a few hanging threads for me, much appreciated. <div><br></div><div>Cheers, Mike</div><div><br></div><div><br></div></div><br><div class="gmail_quote"><div dir="ltr" class="gmail_attr">On Mon, Jun 3, 2024 at 3:48 AM Even Rouault <<a href="mailto:even.rouault@spatialys.com">even.rouault@spatialys.com</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><u></u>
<div>
<p>Michael,</p>
<p><a href="https://github.com/OSGeo/gdal/pull/10108" target="_blank">https://github.com/OSGeo/gdal/pull/10108</a> will fix it.</p>
<p>You can also workaround the issue by overriding the extent of the
source dataset to be exactly -179.995,89.995,180.005,-89.995 with<br>
</p>
<p> "vrt://NETCDF:20240102090000-JPL-L4_GHRSST-SSTfnd-MUR-GLOB-v02.0-fv04.1.nc:analysed_sst?a_srs=EPSG:4326&a_ullr=-179.995,89.995,180.005,-89.995"</p>
<p>The issue was twofolds:</p>
<p>- the dataset uses lon/lat arrays with single precision
floating-point numbers. Hence -179.99 gets actually read a a
read as -179.99000549316406 as a double-precision number, and thus
the maxx - minxx difference was slightly over 360 degrees. I've
added a heuristics to try to round -179.99-single-precision as
-179.99-double-precision<br>
</p>
<p>- which defeated a heuristics in the warper to identify the
center longitude of a dataset registered in a geographic CRS. This
center longitude is just (min_lon + max_lon) / 2 if the dataset
doesn't cover more than 360 degrees of longitude. This value, when
computed, is passed as a hint OGRCoordinateTransformation so that
it can post-correct longitudes to apply a +/- 360 degree offset,
to be in the range of the source dataset. I've relaxed the sanity
check to allow slighly more than 360 degrees</p>
<p>Even<br>
</p>
<p><br>
</p>
<div>Le 01/06/2024 à 08:19, Michael Sumner
via gdal-dev a écrit :<br>
</div>
<blockquote type="cite">
<div dir="ltr">This data source has an odd georeferencing, it's a
36000x17999 raster in the extent -179.995 -89.995 180.005
89.995. <br>
<br>
vrt://NETCDF:/vsicurl/<a href="https://archive.podaac.earthdata.nasa.gov/podaac-ops-cumulus-protected/MUR-JPL-L4-GLOB-v4.1/20240101090000-JPL-L4_GHRSST-SSTfnd-MUR-GLOB-v02.0-fv04.1.nc:analysed_sst?a_srs=EPSG:4326" target="_blank">https://archive.podaac.earthdata.nasa.gov/podaac-ops-cumulus-protected/MUR-JPL-L4-GLOB-v4.1/20240101090000-JPL-L4_GHRSST-SSTfnd-MUR-GLOB-v02.0-fv04.1.nc:analysed_sst?a_srs=EPSG:4326</a><br>
<br>
(Why is it 17999: well there's one missing row - you'd expect
18000 at 0.01 resolution and that manifests as half a cell short
at the top and bottom).
<div><br>
</div>
<div>But also the x registration is half a cell to the east from
what you would expect. <br>
<br>
It's correctly registered in that frame though, if we shift it
west by the half cell I've verified visually that it "looks
wrong" compared to coastlines etc. <br>
<br>
The point of this report is, gdalwarp misses the western
window of data for a request across the antimeridian (there's
no data in the output east of 180). <br>
<br>
gdalwarp $dsn -te 173.00 -19.85 185.00 -15.00 out1.tif -co
COMPRESS=LZW<br>
<br>
If we "fix" the extent to be more longitudinally pleasing, we
get the data east of the antimeridian, the map is properly
filled either side of the antimeridian. <br>
<br>
<br>
vrt://NETCDF:/vsicurl/<a href="https://archive.podaac.earthdata.nasa.gov/podaac-ops-cumulus-protected/MUR-JPL-L4-GLOB-v4.1/20240101090000-JPL-L4_GHRSST-SSTfnd-MUR-GLOB-v02.0-fv04.1.nc:analysed_sst?a_srs=EPSG:4326&a_ullr=-180,89.995,180,-89.995" target="_blank">https://archive.podaac.earthdata.nasa.gov/podaac-ops-cumulus-protected/MUR-JPL-L4-GLOB-v4.1/20240101090000-JPL-L4_GHRSST-SSTfnd-MUR-GLOB-v02.0-fv04.1.nc:analysed_sst?a_srs=EPSG:4326&a_ullr=-180,89.995,180,-89.995</a><br>
<br>
gdalwarp $dsnfix -te 173.00 -19.85 185.00 -15.00 out2.tif -co
COMPRESS=LZW<br>
<br>
<br>
Which is good, but as mentioned above the workaround means we
are now half a cell out. It works fine for projected windows
over the antimeridian in either situation.<br>
<br>
The compressed file sizes reflect the missing data in the
second one, these are 625K and 1.1Mb. <br>
<br>
I'd like the warper to be able to correctly fill the data from
the western window in the original case, we have determined
that we can't feasibly "fix" the extent. </div>
<div><br>
</div>
<div>Or maybe something else I'm missing? I don't *think*
SAMPLE_GRID or SOURCE_EXTRA helps here. I know that I could
craft a meridian-crossing combination in VRT but I'd rather
like this workflow to work as-is. (Chasing down the weird
grid referencing is another project, but it's a really massive
dataset so I'm pretty unconfident of that occurring). <br>
<br>
Cheers, Mike<br clear="all">
<div><br>
</div>
<span class="gmail_signature_prefix">-- </span><br>
<div dir="ltr" class="gmail_signature">
<div dir="ltr">Michael Sumner<br>
Research Software Engineer<br>
Australian Antarctic Division<br>
Hobart, Australia<br>
e-mail: <a href="mailto:mdsumner@gmail.com" target="_blank">mdsumner@gmail.com</a></div>
</div>
</div>
</div>
<br>
<fieldset></fieldset>
<pre>_______________________________________________
gdal-dev mailing list
<a href="mailto:gdal-dev@lists.osgeo.org" target="_blank">gdal-dev@lists.osgeo.org</a>
<a href="https://lists.osgeo.org/mailman/listinfo/gdal-dev" target="_blank">https://lists.osgeo.org/mailman/listinfo/gdal-dev</a>
</pre>
</blockquote>
<pre cols="72">--
<a href="http://www.spatialys.com" target="_blank">http://www.spatialys.com</a>
My software is free, but my time generally not.</pre>
</div>
</blockquote></div><br clear="all"><div><br></div><span class="gmail_signature_prefix">-- </span><br><div dir="ltr" class="gmail_signature"><div dir="ltr">Michael Sumner<br>Research Software Engineer<br>Australian Antarctic Division<br>Hobart, Australia<br>e-mail: <a href="mailto:mdsumner@gmail.com" target="_blank">mdsumner@gmail.com</a></div></div>