Opened 11 years ago

Closed 9 years ago

#4891 closed defect (fixed)

HDF4 scale/offset CoastWatch convention questionable

Reported by: warmerdam Owned by: warmerdam
Priority: normal Milestone: 1.10.2
Component: GDAL_Raster Version: unspecified
Severity: normal Keywords: hdf4
Cc:

Description

The offset and scale reported for:

http://download.osgeo.org/gdal/data/hdf4/REANALYSIS_1999217.hdf

are:

gdalinfo HDF4_SDS:UNKNOWN:"REANALYSIS_1999217.hdf":3 
Driver: HDF4Image/HDF4 Dataset
Files: REANALYSIS_1999217.hdf
Size is 144, 73
Coordinate System is `'
Metadata:
  actual_range=-6.800000191, 78.80000305
  add_offset=277.6499939
  base_date=1999, 1, 1
  Conventions=COARDS
  dataset=NMC Reanalysis
  Day Of Year=217
  description=Data is from NMC initialized reanalysis
(4x/day).  It consists of most variables interpolated to
pressure surfaces from model (sigma) surfaces.
  GRIB_id=54
  GRIB_name=PWAT
  history=created 99/01/04 by Hoop (netCDF2.3)
  least_significant_digit=-1
  level_desc=Entire Atmosphere Considered As a Single Layer
  long_name=4xDaily Precipitable Water for entire atmosphere
  missing_value=32766
  parent_stat=Other
  platform=Model
  precision=2
  references=http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.html
  scale_factor=0.009999999776
  statistic=Individual Obs
  title=4x daily NMC reanalysis (1999)
  units=kg/m^2
  unpacked_valid_range=-50, 150
  valid_range=-32765, -12765
  var_desc=Precipitable Water Content
Corner Coordinates:
...
Band 1 Block=144x73 Type=Int16, ColorInterp=Gray
  NoData Value=32766
  Offset: -2.7764998768064,   Scale:0.009999999776
...

In particular the offset is -2.7765 but it *should* be 277.65. Reviewing the code we seem to be going through:

          // Coastwatch offset and scale.
          if( CSLFetchNameValue( poDS->papszLocalMetadata, "scale_factor" )
              && CSLFetchNameValue( poDS->papszLocalMetadata, "add_offset" ) )
          {
              for( i = 1; i <= poDS->nBands; i++ )
              {
                  HDF4ImageRasterBand *poBand =
                      (HDF4ImageRasterBand *) poDS->GetRasterBand(i);

                  poBand->bHaveScale = poBand->bHaveOffset = TRUE;
                  poBand->dfScale =
                      CPLAtof( CSLFetchNameValue( poDS->papszLocalMetadata,
                                                  "scale_factor" ) );
                  poBand->dfOffset = -1 * poBand->dfScale *
                      CPLAtof( CSLFetchNameValue( poDS->papszLocalMetadata,
                                                  "add_offset" ) );
              }
          }

I believe the reverse of sign and scaling of the offset is inappropriate to some, most or possibly even all products of this style.

Change History (3)

comment:1 by warmerdam, 11 years ago

The CoastWatch logic was added in r8800 by myself, but I can't seem to find any of my CoastWatch data stash that actually has a non-zero add_offset value. Based on this finding I am going to use a more conventional logic for the offset.

comment:2 by warmerdam, 11 years ago

Status: newassigned

Change made in trunk (r25219).

Should consider back porting to 1.9 branch.

comment:3 by Jukka Rahkonen, 9 years ago

Milestone: 1.10.2
Resolution: fixed
Status: assignedclosed

Change has been living in GDAL code for two years without complaint.

Note: See TracTickets for help on using tickets.