Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

netCDF: use the SRS (its geographic part) if found in the file,… #9528

Merged
merged 1 commit into from
Mar 25, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 8 additions & 4 deletions autotest/gdrivers/netcdf.py
Original file line number Diff line number Diff line change
Expand Up @@ -1576,11 +1576,13 @@ def test_netcdf_42():
gdaltest.netcdf_drv.CreateCopy("tmp/netcdf_42.nc", src_ds)

ds = gdal.Open("tmp/netcdf_42.nc")
assert ds.GetMetadata("GEOLOCATION") == {
got_md = ds.GetMetadata("GEOLOCATION")
assert got_md["SRS"].startswith('GEOGCRS["WGS 84",')
del got_md["SRS"]
assert got_md == {
"LINE_OFFSET": "0",
"X_DATASET": 'NETCDF:"tmp/netcdf_42.nc":lon',
"PIXEL_STEP": "1",
"SRS": 'GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]',
"PIXEL_OFFSET": "0",
"X_BAND": "1",
"LINE_STEP": "1",
Expand Down Expand Up @@ -1718,11 +1720,13 @@ def test_netcdf_43():
)

ds = gdal.Open("tmp/netcdf_43.nc")
assert ds.GetMetadata("GEOLOCATION") == {
got_md = ds.GetMetadata("GEOLOCATION")
assert got_md["SRS"].startswith('GEOGCRS["NAD27",')
del got_md["SRS"]
assert got_md == {
"LINE_OFFSET": "0",
"X_DATASET": 'NETCDF:"tmp/netcdf_43.nc":lon',
"PIXEL_STEP": "1",
"SRS": 'GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AXIS["Latitude",NORTH],AXIS["Longitude",EAST],AUTHORITY["EPSG","4326"]]',
"PIXEL_OFFSET": "0",
"X_BAND": "1",
"LINE_STEP": "1",
Expand Down
25 changes: 21 additions & 4 deletions frmts/netcdf/netcdfdataset.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -4266,9 +4266,25 @@ void netCDFDataset::SetProjectionFromVar(
if (bReadSRSOnly)
return;

// Determines the SRS to be used by the geolocation array, if any
std::string osGeolocWKT = SRS_WKT_WGS84_LAT_LONG;
if (!m_oSRS.IsEmpty())
{
OGRSpatialReference oGeogCRS;
oGeogCRS.CopyGeogCSFrom(&m_oSRS);
char *pszWKTTmp = nullptr;
const char *const apszOptions[] = {"FORMAT=WKT2_2019", nullptr};
if (oGeogCRS.exportToWkt(&pszWKTTmp, apszOptions) == OGRERR_NONE)
{
osGeolocWKT = pszWKTTmp;
}
CPLFree(pszWKTTmp);
}

// Process geolocation arrays from CF "coordinates" attribute.
std::string osGeolocXName, osGeolocYName;
if (ProcessCFGeolocation(nGroupId, nVarId, osGeolocXName, osGeolocYName))
if (ProcessCFGeolocation(nGroupId, nVarId, osGeolocWKT, osGeolocXName,
osGeolocYName))
{
bool bCanCancelGT = true;
if ((nVarDimXID != -1) && (nVarDimYID != -1))
Expand Down Expand Up @@ -4311,7 +4327,7 @@ void netCDFDataset::SetProjectionFromVar(
CPLDebug("GDAL_netCDF", "using variables %s and %s for GEOLOCATION",
pszGeolocXFullName, pszGeolocYFullName);

GDALPamDataset::SetMetadataItem("SRS", SRS_WKT_WGS84_LAT_LONG,
GDALPamDataset::SetMetadataItem("SRS", osGeolocWKT.c_str(),
"GEOLOCATION");

CPLString osTMP;
Expand Down Expand Up @@ -4729,6 +4745,7 @@ bool netCDFDataset::ProcessNASAEMITGeoLocation(int nGroupId, int nVarId)
}

int netCDFDataset::ProcessCFGeolocation(int nGroupId, int nVarId,
const std::string &osGeolocWKT,
std::string &osGeolocXNameOut,
std::string &osGeolocYNameOut)
{
Expand Down Expand Up @@ -4801,8 +4818,8 @@ int netCDFDataset::ProcessCFGeolocation(int nGroupId, int nVarId,
"using variables %s and %s for GEOLOCATION",
pszGeolocXFullName, pszGeolocYFullName);

GDALPamDataset::SetMetadataItem(
"SRS", SRS_WKT_WGS84_LAT_LONG, "GEOLOCATION");
GDALPamDataset::SetMetadataItem("SRS", osGeolocWKT.c_str(),
"GEOLOCATION");

CPLString osTMP;
osTMP.Printf("NETCDF:\"%s\":%s", osFilename.c_str(),
Expand Down
1 change: 1 addition & 0 deletions frmts/netcdf/netcdfdataset.h
Original file line number Diff line number Diff line change
Expand Up @@ -502,6 +502,7 @@ class netCDFDataset final : public GDALPamDataset
bool ProcessNASAEMITGeoLocation(int nGroupId, int nVarId);

int ProcessCFGeolocation(int nGroupId, int nVarId,
const std::string &osGeolocWKT,
std::string &osGeolocXNameOut,
std::string &osGeolocYNameOut);
CPLErr Set1DGeolocation(int nGroupId, int nVarId, const char *szDimName);
Expand Down
Loading