ucar.nc2.dataset.conv.HdfEosModisConvention Maven / Gradle / Ivy
The newest version!
/*
* Copyright (c) 1998-2018 John Caron and University Corporation for Atmospheric Research/Unidata
* See LICENSE for license information.
*/
package ucar.nc2.dataset.conv;
import ucar.ma2.Array;
import ucar.ma2.DataType;
import ucar.nc2.*;
import ucar.nc2.constants.*;
import ucar.nc2.dataset.*;
import ucar.nc2.iosp.hdf4.HdfEos;
import ucar.nc2.time.CalendarDate;
import ucar.nc2.util.CancelTask;
import ucar.unidata.geoloc.projection.Sinusoidal;
/**
* HDF4-EOS TERRA MODIS
*
* @author caron
* @since 2/23/13
* @see "https://lpdaac.usgs.gov/products/modis_overview"
*/
public class HdfEosModisConvention extends ucar.nc2.dataset.CoordSysBuilder {
private static final String CRS = "Projection";
private static final String DATA_GROUP = "Data_Fields";
private static final String DIMX_NAME = "XDim";
private static final String DIMY_NAME = "YDim";
private static final String TIME_NAME = "time";
public static boolean isMine(NetcdfFile ncfile) {
if (!ncfile.getFileTypeId().equals("HDF4-EOS"))
return false;
String typeName = ncfile.findAttValueIgnoreCase(null, CF.FEATURE_TYPE, null);
if (typeName == null)
return false;
if (!typeName.equals(FeatureType.GRID.toString()) && !typeName.equals(FeatureType.SWATH.toString()))
return false;
return checkGroup(ncfile.getRootGroup());
}
private static boolean checkGroup(Group g) {
Variable crs = g.findVariableLocal(HdfEos.HDFEOS_CRS);
Group dataG = g.findGroupLocal(DATA_GROUP);
if (crs != null && dataG != null) {
Attribute att = crs.findAttribute(HdfEos.HDFEOS_CRS_Projection);
if (att == null)
return false;
if (!att.getStringValue().equals("GCTP_SNSOID") && !att.getStringValue().equals("GCTP_GEO"))
return false;
return !(dataG.findDimensionLocal(DIMX_NAME) == null || dataG.findDimensionLocal(DIMY_NAME) == null);
}
for (Group ng : g.getGroups()) {
if (checkGroup(ng))
return true;
}
return false;
}
public HdfEosModisConvention() {
this.conventionName = "HDF4-EOS-MODIS";
}
/*
* group: MODIS_Grid_16DAY_250m_500m_VI {
* variables:
* short _HDFEOS_CRS;
* :Projection = "GCTP_SNSOID";
* :UpperLeftPointMtrs = -2.0015109354E7, 1111950.519667; // double
* :LowerRightMtrs = -1.8903158834333E7, 0.0; // double
* :ProjParams = 6371007.181, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0; // double
* :SphereCode = "-1";
*
*
* group: Data_Fields {
* dimensions:
* YDim = 4800;
* XDim = 4800;
* variables:
* short 250m_16_days_NDVI(YDim=4800, XDim=4800);
* ...
*
*/
private boolean addTimeCoord;
public void augmentDataset(NetcdfDataset ds, CancelTask cancelTask) {
addTimeCoord = addTimeCoordinate(ds);
augmentGroup(ds, ds.getRootGroup());
ds.addAttribute(ds.getRootGroup(), new Attribute(CDM.CONVENTIONS, "CF-1.0"));
ds.finish();
}
private boolean addTimeCoordinate(NetcdfDataset ds) {
// add time coordinate
CalendarDate cd = parseFilenameForDate(ds.getLocation());
if (cd == null)
return false;
ds.addAttribute(ds.getRootGroup(), new Attribute("_MODIS_Date", cd.toString()));
// add the time dimension
int nTimesDim = 1;
Dimension newDim = new Dimension(TIME_NAME, nTimesDim);
ds.addDimension(null, newDim);
// add the coordinate variable
String units = "seconds since " + cd;
String desc = "time coordinate";
Array data = Array.makeArray(DataType.DOUBLE, 1, 0.0, 0.0);
CoordinateAxis1D timeCoord = new CoordinateAxis1D(ds, null, TIME_NAME, DataType.DOUBLE, "", units, desc);
timeCoord.setCachedData(data, true);
timeCoord.addAttribute(new Attribute(_Coordinate.AxisType, AxisType.Time.toString()));
ds.addCoordinateAxis(timeCoord);
return true;
}
private CalendarDate parseFilenameForDate(String filename) {
// filename MOD13Q1.A2000065.h11v04.005.2008238031620.hdf
String[] tokes = filename.split("\\.");
if (tokes.length < 2)
return null;
if (tokes[1].length() < 8)
return null;
String want = tokes[1];
String yearS = want.substring(1, 5);
String jdayS = want.substring(5, 8);
try {
int year = Integer.parseInt(yearS);
int jday = Integer.parseInt(jdayS);
return CalendarDate.withDoy(null, year, jday, 0, 0, 0);
} catch (Exception e) {
return null;
}
}
private void augmentGroup(NetcdfDataset ds, Group g) {
Variable crs = g.findVariableLocal(HdfEos.HDFEOS_CRS);
if (crs != null)
augmentGroupWithProjectionInfo(ds, g);
for (Group ng : g.getGroups())
augmentGroup(ds, ng);
}
private void augmentGroupWithProjectionInfo(NetcdfDataset ds, Group g) {
Dimension dimX = null, dimY = null;
Group dataG = g.findGroupLocal(DATA_GROUP);
if (dataG != null) {
dimX = dataG.findDimensionLocal(DIMX_NAME);
dimY = dataG.findDimensionLocal(DIMY_NAME);
}
if (dimX == null || dimY == null)
return;
Variable crs = g.findVariableLocal(HdfEos.HDFEOS_CRS);
Attribute projAtt = crs.findAttribute(HdfEos.HDFEOS_CRS_Projection);
if (projAtt != null) {
Attribute upperLeft = crs.findAttribute(HdfEos.HDFEOS_CRS_UpperLeft);
Attribute lowerRight = crs.findAttribute(HdfEos.HDFEOS_CRS_LowerRight);
Attribute projParams = crs.findAttribute(HdfEos.HDFEOS_CRS_ProjParams);
double minX = upperLeft.getNumericValue(0).doubleValue();
double minY = upperLeft.getNumericValue(1).doubleValue();
double maxX = lowerRight.getNumericValue(0).doubleValue();
double maxY = lowerRight.getNumericValue(1).doubleValue();
boolean hasProjection = false;
String coordinates = null;
ProjectionCT ct;
if (projAtt.getStringValue().equals("GCTP_SNSOID")) {
hasProjection = true;
ct = makeSinusoidalProjection(CRS, projParams);
VariableDS crss = makeCoordinateTransformVariable(ds, ct);
crss.addAttribute(new Attribute(_Coordinate.AxisTypes, "GeoX GeoY"));
ds.addVariable(dataG, crss);
ds.addCoordinateAxis(makeCoordAxis(ds, dataG, DIMX_NAME, dimX.getLength(), minX, maxX, true));
ds.addCoordinateAxis(makeCoordAxis(ds, dataG, DIMY_NAME, dimY.getLength(), minY, maxY, false));
coordinates = addTimeCoord ? TIME_NAME + " " + DIMX_NAME + " " + DIMY_NAME : DIMX_NAME + " " + DIMY_NAME;
} else if (projAtt.getStringValue().equals("GCTP_GEO")) {
ds.addCoordinateAxis(makeLatLonCoordAxis(ds, dataG, dimX.getLength(), minX * 1e-6, maxX * 1e-6, true));
ds.addCoordinateAxis(makeLatLonCoordAxis(ds, dataG, dimY.getLength(), minY * 1e-6, maxY * 1e-6, false));
coordinates = addTimeCoord ? TIME_NAME + " Lat Lon" : "Lat Lon";
}
for (Variable v : dataG.getVariables()) {
if (v.getRank() != 2)
continue;
if (!v.getDimension(0).equals(dimY))
continue;
if (!v.getDimension(1).equals(dimX))
continue;
if (coordinates != null)
v.addAttribute(new Attribute(CF.COORDINATES, coordinates));
if (hasProjection)
v.addAttribute(new Attribute(CF.GRID_MAPPING, CRS));
}
}
}
/*
* The UpperLeftPointMtrs is in projection coordinates, and identifies the very upper left corner of the upper left
* pixel of the image data
* • The LowerRightMtrs identifies the very lower right corner of the lower right pixel of the image data. These
* projection coordinates are the only metadata that accurately reflect the extreme corners of the gridded image
*/
private CoordinateAxis makeCoordAxis(NetcdfDataset ds, Group g, String name, int n, double start, double end,
boolean isX) {
CoordinateAxis v =
new CoordinateAxis1D(ds, g, name, DataType.DOUBLE, name, "km", isX ? "x coordinate" : "y coordinate");
double incr = (end - start) / n;
v.setValues(n, start * .001, incr * .001); // km
v.addAttribute(new Attribute(_Coordinate.AxisType, isX ? AxisType.GeoX.toString() : AxisType.GeoY.toString()));
return v;
}
/*
* group: MOD_Grid_MOD17A3 {
* variables:
* short _HDFEOS_CRS;
* :Projection = "GCTP_GEO";
* :UpperLeftPointMtrs = -1.8E8, 9.0E7; // double
* :LowerRightMtrs = 1.8E8, -9.0E7; // double
*/
private CoordinateAxis makeLatLonCoordAxis(NetcdfDataset ds, Group g, int n, double start, double end,
boolean isLon) {
String name = isLon ? AxisType.Lon.toString() : AxisType.Lat.toString();
String dimName = isLon ? DIMX_NAME : DIMY_NAME;
CoordinateAxis v =
new CoordinateAxis1D(ds, g, name, DataType.DOUBLE, dimName, isLon ? "degrees_east" : "degrees_north", null);
double incr = (end - start) / n;
v.setValues(n, start, incr);
v.addAttribute(new Attribute(_Coordinate.AxisType, name));
return v;
}
/*
* 1 5 7 8
* 16 PGSd_SNSOID Sphere CentMer FE FN
*
*/
private ProjectionCT makeSinusoidalProjection(String name, Attribute projParams) {
double radius = projParams.getNumericValue(0).doubleValue();
double centMer = projParams.getNumericValue(4).doubleValue();
double falseEast = projParams.getNumericValue(6).doubleValue();
double falseNorth = projParams.getNumericValue(7).doubleValue();
Sinusoidal proj = new Sinusoidal(centMer, falseEast * .001, falseNorth * .001, radius * .001);
return new ProjectionCT(name, "FGDC", proj);
}
}
/*
*
* Can we extract projection info ??
*
* // external xml
*
* Terra
*
* MODIS
*
* MODIS
*
*
*
*
* // ODL
* UpperLeftPointMtrs=(-20015109.354000,1111950.519667)
* LowerRightMtrs=(-18903158.834333,-0.000000)
* Projection=GCTP_SNSOID
* ProjParams=(6371007.181000,0,0,0,0,0,0,0,0,0,0,0,0)
* SphereCode=-1
*
* // from http://datamirror.csdb.cn/modis/resource/doc/HDFEOS5ref.pdf
*
* 1 5 7 8
* 16 PGSd_SNSOID Sphere CentMer FE FN
*
* CentMer Longitude of the central meridian
* OriginLat Latitude of the projection origin
* FE False easting in the same units as the semi-major axis
* FN False northing in the same units as the semi-major axis
*
* according to snyder (p 243), sinusoidal projection parameterized only by CentMer and R
* its a very simple form for sphere.
*
* question as usual is what are the x,y coordinates? Assuming even spacing from the cornerss is the only thing you can
* do:
*
* UpperLeftPointMtrs=(-20015109.354000, 1111950.519667)
* LowerRightMtrs=(-18903158.834333,-0.000000)
*
*
*/
© 2015 - 2025 Weber Informatics LLC | Privacy Policy