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

[5.5.3]: Fix GRIB IOSP time coords from intervals, for 2-D time datasets #1248

Merged
Show file tree
Hide file tree
Changes from 2 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
49 changes: 5 additions & 44 deletions grib/src/main/java/ucar/nc2/grib/collection/GribIosp.java
Original file line number Diff line number Diff line change
Expand Up @@ -776,58 +776,18 @@ private Array makeLazyTime2Darray(Variable coord, Time2Dinfo info) {
break;

case intv:
for (int runIdx = 0; runIdx < nruns; runIdx++) {
CoordinateTimeIntv timeIntv = (CoordinateTimeIntv) time2D.getTimeCoordinate(runIdx);
int timeIdx = 0;
for (TimeCoordIntvValue tinv : timeIntv.getTimeIntervals()) {
data[runIdx * ntimes + timeIdx] = timeUnit.getValue() * tinv.getCoordValue() + time2D.getOffset(runIdx);
timeIdx++;
}
}
break;

case intvU:
// data[nruns*ntimes]
int dataIndex = 0;
double[] currentBounds = new double[2];
double[] prevBounds = new double[2];
count = 0; // dataIndex
for (int runIdx = 0; runIdx < nruns; runIdx++) {
CoordinateTimeIntv timeIntv = (CoordinateTimeIntv) time2D.getTimeCoordinate(runIdx);
int runOffsetIndex = runIdx * ntimes;
int timeUnitValue = timeUnit.getValue();
int time2D_Offset = time2D.getOffset(runIdx);
for (TimeCoordIntvValue tinv : timeIntv.getTimeIntervals()) {
currentBounds[0] = timeUnitValue * tinv.getBounds1() + time2D_Offset;
currentBounds[1] = timeUnitValue * tinv.getBounds2() + time2D_Offset;
// Use end-point of current interval as initial guess for current time coordinate value
data[dataIndex] = currentBounds[1];
if (dataIndex >= runOffsetIndex + 1) {
// Check that time coordinate values are increasing in a strictly-monotonic manner
// (as required by CF conventions). If not strictly-monotonic ...
if (data[dataIndex] <= data[dataIndex - 1]) {
if (dataIndex >= runOffsetIndex + 2) {
if (data[dataIndex - 2] <= currentBounds[0]) {
// Change previous time coordinate value to mid-point between
// current time interval start and end values.
data[dataIndex - 1] = (currentBounds[1] - currentBounds[0]) / 2.0 + currentBounds[0];
} else {
// Or change previous time coordinate value to mid-point between
// current time interval end value and the time coord value from two steps back.
data[dataIndex - 1] = (currentBounds[1] - data[dataIndex - 2]) / 2.0 + data[dataIndex - 2];
}
} else {
data[dataIndex - 1] = (prevBounds[1] - prevBounds[0]) / 2.0 + prevBounds[0];
}
}
}
prevBounds[0] = currentBounds[0];
prevBounds[1] = currentBounds[1];
dataIndex++;
}
assert (count == runOffsetIndex);
count = GribTimeCoordIntervalUtils.generateTimeCoordValuesFromTimeCoordIntervals(timeIntv.getTimeIntervals(),
data, null, runOffsetIndex, timeUnitValue, time2D_Offset);
}
// The above assumes that sets of intervals are always sorted by
// end point then starting point. That sorting scheme is implemented in
// ucar.nc2.grib.coord.TimeCoordIntvValue.compareTo(o) -- 21 Dec 2022.
break;

case is1Dtime:
Expand Down Expand Up @@ -881,6 +841,7 @@ private Array makeLazyTime2Darray(Variable coord, Time2Dinfo info) {
return Array.factory(DataType.DOUBLE, coord.getShape(), data);
}


private void makeTimeCoordinate1D(NetcdfFile ncfile, Group g, CoordinateTime coordTime) { // }, CoordinateRuntime
// runtime) {
int ntimes = coordTime.getSize();
Expand Down
21 changes: 10 additions & 11 deletions grib/src/main/java/ucar/nc2/grib/collection/GribIospBuilder.java
Original file line number Diff line number Diff line change
Expand Up @@ -676,26 +676,25 @@ private String makeTimeOffsetOrthogonal(Group.Builder g, CoordinateTime2D time2D
v.addAttribute(new Attribute(CDM.UDUNITS, time2D.getTimeUdUnit()));
v.addAttribute(new Attribute(_Coordinate.AxisType, AxisType.TimeOffset.toString()));

double[] midpoints = new double[n];
double[] data = new double[n];
double[] bounds = null;
int timeUnitValue = time2D.getTimeUnit().getValue();

if (time2D.isTimeInterval()) {
bounds = new double[2 * n];
int count = 0;
int countb = 0;
for (Object offset : offsets) {
TimeCoordIntvValue tinv = (TimeCoordIntvValue) offset;
midpoints[count++] = (tinv.getBounds1() + tinv.getBounds2()) / 2.0;
bounds[countb++] = tinv.getBounds1();
bounds[countb++] = tinv.getBounds2();
}
List<TimeCoordIntvValue> intervals = (List<TimeCoordIntvValue>) offsets;
// TODO Should this be using 'timeUnitValue'? Wasn't in previous incarnation.
int count = GribTimeCoordIntervalUtils.generateTimeCoordValuesFromTimeCoordIntervals(intervals, data, bounds, 0,
timeUnitValue, 0);
assert (count == n);
} else {
int count = 0;
for (Object val : offsets) {
Integer off = (Integer) val;
midpoints[count++] = off; // int ??
data[count++] = timeUnitValue * off; // int ??
}
}
v.setCachedData(Array.factory(DataType.DOUBLE, new int[] {n}, midpoints), false);
v.setCachedData(Array.factory(DataType.DOUBLE, new int[] {n}, data), false);

if (time2D.isTimeInterval()) {
String boundsName = toName + "_bounds";
Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1,75 @@
package ucar.nc2.grib.collection;

import ucar.nc2.grib.coord.TimeCoordIntvValue;

import java.util.List;

public class GribTimeCoordIntervalUtils {
ethanrd marked this conversation as resolved.
Show resolved Hide resolved
/**
* Default location in a time coordinate interval to use for the corresponding time coordinate value
* when generating a time coordinate variable.
*/
enum TimeCoordValueDefaultLocationInIntervals {
endPoint, midPoint, startingPoint
}

ethanrd marked this conversation as resolved.
Show resolved Hide resolved
/**
* Generate values for 1-D time coordinate variable given an ordered set of time coordinate intervals
* (supports datasets with 2D run/forecast time),
* ensures that resulting values are increasing in a strictly-monotonic manner.
*
* Assumes that sets of intervals are always sorted by
* end point then starting point. That sorting scheme is implemented in
* ucar.nc2.grib.coord.TimeCoordIntvValue.compareTo(o) -- 21 Dec 2022.
*
* @param timeIntervals given list of time coord intervals from which to calculate time coordinate values
* @param timeCoordValues array for storing generated time coord values, 1D array that may represent a 2D data array
* (run/forecast)
* @param bounds array for storing interval bounds, may be null if bounds are loaded elsewhere.
* @param runOffsetIndex
* @param timeUnitValue time unit
* @param time2D_Offset
* @return the index into the timeCoordValues array following the last value written
*/
public static int generateTimeCoordValuesFromTimeCoordIntervals(List<TimeCoordIntvValue> timeIntervals,
double[] timeCoordValues, double[] bounds, int runOffsetIndex, int timeUnitValue, int time2D_Offset) {
int dataIndex = runOffsetIndex;
int boundsIndex = runOffsetIndex * 2;
double[] currentBounds = new double[2];
double[] prevBounds = new double[2];
for (TimeCoordIntvValue tinv : timeIntervals) {
currentBounds[0] = timeUnitValue * tinv.getBounds1() + time2D_Offset;
currentBounds[1] = timeUnitValue * tinv.getBounds2() + time2D_Offset;
// Use end-point of current interval as initial guess for current time coordinate value
timeCoordValues[dataIndex] = currentBounds[1];
if (dataIndex >= runOffsetIndex + 1) {
// Check that time coordinate values are increasing in a strictly-monotonic manner
// (as required by CF conventions). If not strictly-monotonic ...
if (timeCoordValues[dataIndex] <= timeCoordValues[dataIndex - 1]) {
if (dataIndex >= runOffsetIndex + 2) {
if (timeCoordValues[dataIndex - 2] <= currentBounds[0]) {
// Change previous time coordinate value to mid-point between
// current time interval start and end values.
timeCoordValues[dataIndex - 1] = (currentBounds[1] - currentBounds[0]) / 2.0 + currentBounds[0];
} else {
// Or change previous time coordinate value to mid-point between
// current time interval end value and the time coord value from two steps back.
timeCoordValues[dataIndex - 1] =
(currentBounds[1] - timeCoordValues[dataIndex - 2]) / 2.0 + timeCoordValues[dataIndex - 2];
}
} else {
timeCoordValues[dataIndex - 1] = (prevBounds[1] - prevBounds[0]) / 2.0 + prevBounds[0];
}
}
}
if (bounds != null) {
bounds[boundsIndex++] = currentBounds[0];
bounds[boundsIndex++] = currentBounds[1];
}
prevBounds[0] = currentBounds[0];
prevBounds[1] = currentBounds[1];
dataIndex++;
}
return dataIndex;
}
}
Binary file not shown.
Original file line number Diff line number Diff line change
Expand Up @@ -94,4 +94,66 @@ private void checkTimeCoordinateVariable1D_IsStrictlyMonotonicallyIncreasing(int
prevValue = currentValue;
}
}

@Test
public void testTimeCoord2D_fromIntervals_isStrictlyMonotonicallyIncreasing_Dataset1() throws IOException {
String testfile = "../grib/src/test/data/GFS_Global_onedeg_20220627.TotalPrecip.Out24hrs.grib2";
String varName = "Total_precipitation_surface_Mixed_intervals_Accumulation";

checkTimeCoord2D_fromIntervals_isStrictlyMonotonicallyIncreasing(testfile, varName);
}

private void checkTimeCoord2D_fromIntervals_isStrictlyMonotonicallyIncreasing(String testfile, String varName)
throws IOException {
try (NetcdfFile nc = NetcdfFiles.open(testfile)) {
Variable dataVar = nc.findVariable(varName);
Assert.assertNotNull(dataVar);

Variable reftimeVar = nc.findVariable("reftime");
Array reftimeValues = reftimeVar.read();
tdrwenski marked this conversation as resolved.
Show resolved Hide resolved

Variable timeOffsetVar = nc.findVariable("timeOffset");
Array timeOffsetValues = timeOffsetVar.read();

Variable timeOffsetBoundsVar = nc.findVariable("timeOffset_bounds");
Array timeOffsetBoundsValues = timeOffsetBoundsVar.read();

Variable timeVar = nc.findVariable("time");
Array timeValues = timeVar.read();

Variable timeBoundsVar = nc.findVariable("time_bounds");
Array timeBoundsValues = timeBoundsVar.read();

Dimension timeDim = null;
for (Dimension dim : dataVar.getDimensions()) {
if (dim.getShortName().startsWith("time")) {
timeDim = dim;
break;
}
}
Variable timeCoordVar = nc.findVariable(timeDim.getShortName());
Assert.assertNotNull(timeCoordVar);

Attribute att = timeCoordVar.findAttribute("bounds");
Assert.assertNotNull(att);
Assert.assertEquals(timeDim.getShortName() + "_bounds", att.getStringValue());

Array timeCoordValues = timeCoordVar.read();
checkTimeCoordVariable2D_IsStrictlyMonotonicallyIncreasing(timeDim.getLength(), timeCoordValues);
}
}

private void checkTimeCoordVariable2D_IsStrictlyMonotonicallyIncreasing(int timeDimLength, Array timeCoordValues) {
double currentValue = timeCoordValues.getDouble(0);
double prevValue = currentValue;
StringBuilder valuesSoFar = new StringBuilder();
valuesSoFar.append(currentValue);
for (int i = 1; i < timeDimLength; i++) {
currentValue = timeCoordValues.getDouble(i);
valuesSoFar.append(", ").append(currentValue);
Assert.assertTrue("Not increasing in a strictly-monotonic manner: [" + valuesSoFar + "]",
currentValue > prevValue);
prevValue = currentValue;
}
}
}