TemporalShiftScale: unexpected output time steps
Hello all, and happy new year!
I'm reporting a bug in the Temporal Shift Scale filter, as well as presenting a solution. The bug is clearly in vtk, so I will try to merge the solution into the vtk source. I copied this as issue 16955 in the vtk gitlab.
Versions
Platform: Linux x86_64
I see this bug in:
- ParaView 4.0.1
- ParaView 5.0.1
- ParaView 5.2.0
I also see that in the latest vtk source code (7.1), vtkTemporalShiftScale.cxx is pretty much the same as in the buggy ParaView versions, so I don't think somebody else already handled this.
Bug
foamCase_vtkTemporalShiftScaleBug.tar.gz
The attached example openfoam case contains 10 time steps:
- [0.0, 0.001, 0.002, 0.003, 0.004, 0.005, 0.006, 0.007, 0.008, 0.009]
I open the case in ParaView via the foam.foam file, and apply the TemporalShiftScale filter with the following settings:
- PreShift = 0
- PostShift = 0
- Scale = 1
- Periodic = on
- PeriodicEndCorrection = off
- MaximumNumerOfPeriods = 2
I then expect the following output time-steps:
- [0.000, 0.001, 0.002, 0.003, 0.004, 0.005, 0.006, 0.007, 0.008, 0.009, 0.010, 0.011, 0.012, 0.013, 0.014,0.015, 0.016, 0.017, 0.018, 0.019]
But instead I get the following output:
- [0.000, 0.001, 0.002, 0.003, 0.004, 0.005, 0.006, 0.007, 0.008, 0.009, 0.010, 0.011, 0.012, 0.013, 0.014, 0.015, 0.016, 0.017, 0.018]
When I look at each time-step in the 2nd cycle I see that they don't even correspond how I expect them:
- 0.010 -> 0.001 (expected 0.000)
- 0.011 -> 0.002 (expected 0.001)
- ...
- 0.017 -> 0.008 (expected 0.007)
- 0.018 -> 0.000 (expected 0.008)
Source of bug
In vtkTemporalShiftScale.cxx::RequestInformation() I see that the periodic range is calculated as:
PeriodicRange[0] = OutRange[0] = ForwardConvert(InRange[0]) //= 0.000 in this case
PeriodicRange[1] = OutRange[1] = ForwardConvert(InRange[1]) //= 0.009 in this case
And everywhere the range is needed it is calculated as:
range = PeriodicRange[1] - PeriodicRange[0] //= 0.009 in this case
I expect this range to be 0.010 though.
This is clearly the source of the bug:
- when cycles are repeated, the periodic range is added to the 1st cycle. Thus the following time-steps are actually written: [0.000, 0.001, 0.002, 0.003, 0.004, 0.005, 0.006, 0.007, 0.008, 0.009, 0.009, 0.010, 0.011, 0.012, 0.013, 0.014, 0.015, 0.016, 0.017, 0.018]
- then there is probably some external function throwing away the duplicate 0.009 time step (eg. needed to handle duplicate time-dirs in an openfoam case)
- in the downstream to upstream function vtkTemporalShiftScale.cxx::RequestUpdateExtent(), the modulo-range operator (implemented via floor division and subtraction of that multiple of the range) is applied only to downstream time-steps greater than the range: thus 0.009 refers to 0.009, but 0.018 refers to 0.000
Solution
I implemented the changes proposed below in the attached version of vtkTemporalShiftScale.cxx: vtkTemporalShiftScale.cxx
I will try to follow the official route to get the bugfix merged into the vtk source. (I'm actually just writing this because I'm waiting for the compilation)
deltaT
The solution requires adding the last deltaT to the periodic range:
- so 0.009 becomes 0.009+0.001=0.010 (in case of PeriodicEndCorrection=off)
OutRange[1] also needs to be corrected by subtracting the same deltaT:
- so 0.018 becomes 2*0.010 - 0.001 = 0.019 (both in case of PeriodicEndCorrection=off and PeriodicEndCorrection=on)
Calculation of deltaT:
- deltaT can be determined directly for PeriodicEndCorrection=on by taking the difference of the last two timesteps.
- deltaT can only be guessed for PeriodicEndCorrection=off by taking the average of the previous timesteps: 0.009/(numTimes-1)=0.009/9=0.001 The documentation should warn the user that for non-uniform deltaTs the last deltaT might not be correct, and that he should supply a cyclic copy of the first timestep and set PeriodicEndCorrection=on.
vtkTemporalShiftScale.cxx::RequestUpdateExtent()
The periodic correction to get the appropriate input time, needs to be applied to all times greater or equals to the range end. So:
if (ttime>this->PeriodicRange[1])
{
....
}
becomes:
if (ttime>=this->PeriodicRange[1])
{
....
}
This means that in the example case 0.010 needs to refer to 0.000, and not to 0.010. 0.010 of course doesn't exists in the input, and is subsequently clipped to 0.009, which is not what we want.
FYI: redundant code
vtkTemporalShiftScale.cxx::RequestInformation() has the following redundant lines of code:
else if (this->PeriodicEndCorrection)
{
outTimes[i] = outTimes[o] + m*range;
}
else if (!this->PeriodicEndCorrection)
{
outTimes[i] = outTimes[o] + m*range;
}
I don't know what the author's intent was, so I'm leaving it alone.