Commit 70e6445f authored by Bryan Cazabonne's avatar Bryan Cazabonne

Corrected latitude errors.

parent d2151d93
......@@ -247,7 +247,12 @@ public class EstimatedViennaOneModel implements DiscreteTroposphericModel {
}
// Compute hydrostatique coefficient c
final double coef = ((dofyear - 28) / 365) * 2 * FastMath.PI + psi;
double t0 = 28;
if (latitude < 0) {
// southern hemisphere: t0 = 28 + an integer half of year
t0 += 183;
}
final double coef = ((dofyear - t0) / 365) * 2 * FastMath.PI + psi;
final double ch = c0h + ((FastMath.cos(coef) + 1) * (c11h / 2) + c10h) * (1 - FastMath.cos(latitude));
// General constants | Wet part
......@@ -315,7 +320,12 @@ public class EstimatedViennaOneModel implements DiscreteTroposphericModel {
}
// Compute hydrostatique coefficient c
final T coef = psi.add(((dofyear - 28) / 365) * 2 * FastMath.PI);
double t0 = 28;
if (latitude < 0) {
// southern hemisphere: t0 = 28 + an integer half of year
t0 += 183;
}
final T coef = psi.add(((dofyear - t0) / 365) * 2 * FastMath.PI);
final T ch = c11h.divide(2.0).multiply(FastMath.cos(coef).add(1.0)).add(c10h).multiply(1 - FastMath.cos(latitude)).add(c0h);
// General constants | Wet part
......@@ -418,13 +428,14 @@ public class EstimatedViennaOneModel implements DiscreteTroposphericModel {
* @return the height correction, in m
*/
private double computeHeightCorrection(final double elevation, final double height) {
final double fixedHeight = FastMath.max(0.0, height);
final double sinE = FastMath.sin(elevation);
// Ref: Eq. 4
final double function = computeFunction(2.53e-5, 5.49e-3, 1.14e-3, elevation);
// Ref: Eq. 6
final double dmdh = (1 / sinE) - function;
// Ref: Eq. 7
final double correction = dmdh * (height / 1000);
final double correction = dmdh * (fixedHeight / 1000);
return correction;
}
......@@ -444,13 +455,14 @@ public class EstimatedViennaOneModel implements DiscreteTroposphericModel {
*/
private <T extends RealFieldElement<T>> T computeHeightCorrection(final T elevation, final T height, final Field<T> field) {
final T zero = field.getZero();
final T fixedHeight = FastMath.max(zero, height);
final T sinE = FastMath.sin(elevation);
// Ref: Eq. 4
final T function = computeFunction(zero.add(2.53e-5), zero.add(5.49e-3), zero.add(1.14e-3), elevation);
// Ref: Eq. 6
final T dmdh = sinE.reciprocal().subtract(function);
// Ref: Eq. 7
final T correction = dmdh.multiply(height.divide(1000.0));
final T correction = dmdh.multiply(fixedHeight.divide(1000.0));
return correction;
}
......
......@@ -98,7 +98,12 @@ public class GlobalMappingFunctionModel implements MappingFunction {
psi = FastMath.PI;
}
final double coef = ((dofyear + 1 - 28) / 365.25) * 2 * FastMath.PI + psi;
double t0 = 28;
if (latitude < 0) {
// southern hemisphere: t0 = 28 + an integer half of year
t0 += 183;
}
final double coef = ((dofyear + 1 - t0) / 365.25) * 2 * FastMath.PI + psi;
final double ch = c0h + ((FastMath.cos(coef) + 1) * (c11h / 2.0) + c10h) * (1.0 - FastMath.cos(latitude));
// bw and cw constants (Boehm, J et al, 2006) | WET PART
......@@ -180,7 +185,12 @@ public class GlobalMappingFunctionModel implements MappingFunction {
psi = zero.add(FastMath.PI);
}
final T coef = psi.add(((dofyear + 1 - 28) / 365.25) * 2 * FastMath.PI);
double t0 = 28;
if (latitude < 0) {
// southern hemisphere: t0 = 28 + an integer half of year
t0 += 183;
}
final T coef = psi.add(((dofyear + 1 - t0) / 365.25) * 2 * FastMath.PI);
final T ch = c11h.divide(2.0).multiply(FastMath.cos(coef).add(1.0)).add(c10h).multiply(1 - FastMath.cos(latitude)).add(c0h);
// bw and cw constants (Boehm, J et al, 2006) | WET PART
......@@ -291,13 +301,14 @@ public class GlobalMappingFunctionModel implements MappingFunction {
* @return the height correction, in m
*/
private double computeHeightCorrection(final double elevation, final double height) {
final double fixedHeight = FastMath.max(0.0, height);
final double sinE = FastMath.sin(elevation);
// Ref: Eq. 4
final double function = computeFunction(2.53e-5, 5.49e-3, 1.14e-3, elevation);
// Ref: Eq. 6
final double dmdh = (1 / sinE) - function;
// Ref: Eq. 7
final double correction = dmdh * (height / 1000.0);
final double correction = dmdh * (fixedHeight / 1000.0);
return correction;
}
......@@ -317,13 +328,14 @@ public class GlobalMappingFunctionModel implements MappingFunction {
*/
private <T extends RealFieldElement<T>> T computeHeightCorrection(final T elevation, final T height, final Field<T> field) {
final T zero = field.getZero();
final T fixedHeight = FastMath.max(zero, height);
final T sinE = FastMath.sin(elevation);
// Ref: Eq. 4
final T function = computeFunction(zero.add(2.53e-5), zero.add(5.49e-3), zero.add(1.14e-3), elevation);
// Ref: Eq. 6
final T dmdh = sinE.reciprocal().subtract(function);
// Ref: Eq. 7
final T correction = dmdh.multiply(height.divide(1000.0));
final T correction = dmdh.multiply(fixedHeight.divide(1000.0));
return correction;
}
......
......@@ -158,11 +158,20 @@ public class NiellMappingFunctionModel implements MappingFunction {
final int dofyear = dtc.getDate().getDayOfYear();
// Temporal factor
final double coef = 2 * FastMath.PI * ((dofyear - 28) / 365.25);
double t0 = 28;
if (latitude < 0) {
// southern hemisphere: t0 = 28 + an integer half of year
t0 += 183;
}
final double coef = 2 * FastMath.PI * ((dofyear - t0) / 365.25);
final double cosCoef = FastMath.cos(coef);
// Compute ah, bh and ch Eq. 5
final double absLatidude = FastMath.abs(latitude);
double absLatidude = FastMath.abs(latitude);
// there are no data in the model of latitudes lower than 15°
absLatidude = FastMath.max(FastMath.toRadians(15.0), absLatidude);
// there are no data in the model for latitudes greater than 75°
absLatidude = FastMath.min(FastMath.toRadians(75.0), absLatidude);
final double ah = ahAverageFunction.value(absLatidude) - ahAmplitudeFunction.value(absLatidude) * cosCoef;
final double bh = bhAverageFunction.value(absLatidude) - bhAmplitudeFunction.value(absLatidude) * cosCoef;
final double ch = chAverageFunction.value(absLatidude) - chAmplitudeFunction.value(absLatidude) * cosCoef;
......@@ -193,11 +202,20 @@ public class NiellMappingFunctionModel implements MappingFunction {
final int dofyear = dtc.getDate().getDayOfYear();
// Temporal factor
final T coef = zero.add(2 * FastMath.PI * ((dofyear - 28) / 365.25));
double t0 = 28;
if (latitude < 0) {
// southern hemisphere: t0 = 28 + an integer half of year
t0 += 183;
}
final T coef = zero.add(2 * FastMath.PI * ((dofyear - t0) / 365.25));
final T cosCoef = FastMath.cos(coef);
// Compute ah, bh and ch Eq. 5
final double absLatidude = FastMath.abs(latitude);
double absLatidude = FastMath.abs(latitude);
// there are no data in the model of latitudes lower than 15°
absLatidude = FastMath.max(FastMath.toRadians(15.0), absLatidude);
// there are no data in the model for latitudes greater than 75°
absLatidude = FastMath.min(FastMath.toRadians(75.0), absLatidude);
final T ah = cosCoef.multiply(ahAmplitudeFunction.value(absLatidude)).negate().add(ahAverageFunction.value(absLatidude));
final T bh = cosCoef.multiply(bhAmplitudeFunction.value(absLatidude)).negate().add(bhAverageFunction.value(absLatidude));
final T ch = cosCoef.multiply(chAmplitudeFunction.value(absLatidude)).negate().add(chAverageFunction.value(absLatidude));
......@@ -269,13 +287,14 @@ public class NiellMappingFunctionModel implements MappingFunction {
* @return the height correction, in m
*/
private double computeHeightCorrection(final double elevation, final double height) {
final double fixedHeight = FastMath.max(0.0, height);
final double sinE = FastMath.sin(elevation);
// Ref: Eq. 4
final double function = computeFunction(2.53e-5, 5.49e-3, 1.14e-3, elevation);
// Ref: Eq. 6
final double dmdh = (1 / sinE) - function;
// Ref: Eq. 7
final double correction = dmdh * (height / 1000.0);
final double correction = dmdh * (fixedHeight / 1000.0);
return correction;
}
......@@ -289,13 +308,14 @@ public class NiellMappingFunctionModel implements MappingFunction {
*/
private <T extends RealFieldElement<T>> T computeHeightCorrection(final T elevation, final T height, final Field<T> field) {
final T zero = field.getZero();
final T fixedHeight = FastMath.max(zero, height);
final T sinE = FastMath.sin(elevation);
// Ref: Eq. 4
final T function = computeFunction(zero.add(2.53e-5), zero.add(5.49e-3), zero.add(1.14e-3), elevation);
// Ref: Eq. 6
final T dmdh = sinE.reciprocal().subtract(function);
// Ref: Eq. 7
final T correction = dmdh.multiply(height.divide(1000.0));
final T correction = dmdh.multiply(fixedHeight.divide(1000.0));
return correction;
}
......
......@@ -138,8 +138,14 @@ public class ViennaOneModel implements DiscreteTroposphericModel {
psi = FastMath.PI;
}
// Temporal factor
double t0 = 28;
if (latitude < 0) {
// southern hemisphere: t0 = 28 + an integer half of year
t0 += 183;
}
// Compute hydrostatique coefficient c
final double coef = ((dofyear - 28) / 365) * 2 * FastMath.PI + psi;
final double coef = ((dofyear - t0) / 365) * 2 * FastMath.PI + psi;
final double ch = c0h + ((FastMath.cos(coef) + 1) * (c11h / 2) + c10h) * (1 - FastMath.cos(latitude));
// General constants | Wet part
......@@ -187,7 +193,13 @@ public class ViennaOneModel implements DiscreteTroposphericModel {
}
// Compute hydrostatique coefficient c
final T coef = psi.add(((dofyear - 28) / 365) * 2 * FastMath.PI);
// Temporal factor
double t0 = 28;
if (latitude < 0) {
// southern hemisphere: t0 = 28 + an integer half of year
t0 += 183;
}
final T coef = psi.add(((dofyear - t0) / 365) * 2 * FastMath.PI);
final T ch = c11h.divide(2.0).multiply(FastMath.cos(coef).add(1.0)).add(c10h).multiply(1 - FastMath.cos(latitude)).add(c0h);
// General constants | Wet part
......@@ -263,13 +275,14 @@ public class ViennaOneModel implements DiscreteTroposphericModel {
* @return the height correction, in m
*/
private double computeHeightCorrection(final double elevation, final double height) {
final double fixedHeight = FastMath.max(0.0, height);
final double sinE = FastMath.sin(elevation);
// Ref: Eq. 4
final double function = computeFunction(2.53e-5, 5.49e-3, 1.14e-3, elevation);
// Ref: Eq. 6
final double dmdh = (1 / sinE) - function;
// Ref: Eq. 7
final double correction = dmdh * (height / 1000);
final double correction = dmdh * (fixedHeight / 1000);
return correction;
}
......@@ -289,13 +302,14 @@ public class ViennaOneModel implements DiscreteTroposphericModel {
*/
private <T extends RealFieldElement<T>> T computeHeightCorrection(final T elevation, final T height, final Field<T> field) {
final T zero = field.getZero();
final T fixedHeight = FastMath.max(zero, height);
final T sinE = FastMath.sin(elevation);
// Ref: Eq. 4
final T function = computeFunction(zero.add(2.53e-5), zero.add(5.49e-3), zero.add(1.14e-3), elevation);
// Ref: Eq. 6
final T dmdh = sinE.reciprocal().subtract(function);
// Ref: Eq. 7
final T correction = dmdh.multiply(height.divide(1000.0));
final T correction = dmdh.multiply(fixedHeight.divide(1000.0));
return correction;
}
......
......@@ -67,7 +67,7 @@ public class EstimatedViennaOneModelTest {
@Test
public void testZHDStartParameterDerivative() {
doTestParametersDerivatives(EstimatedViennaOneModel.START_HYDROSTATIC_ZENITH_DELAY, 1.0e-16);
doTestParametersDerivatives(EstimatedViennaOneModel.START_HYDROSTATIC_ZENITH_DELAY, 1.0e-14);
}
@Test
......@@ -87,7 +87,7 @@ public class EstimatedViennaOneModelTest {
@Test
public void testAHStartParameterDerivative() {
doTestParametersDerivatives(EstimatedViennaOneModel.START_AH_COEFFICIENT, 7.7e-13);
doTestParametersDerivatives(EstimatedViennaOneModel.START_AH_COEFFICIENT, 3.0e-12);
}
@Test
......@@ -373,7 +373,7 @@ public class EstimatedViennaOneModelTest {
}
for (int i = 0; i < 6; i++) {
Assert.assertEquals(compDeriv[i + 1], refDeriv[0][i], 1.8e-11);
Assert.assertEquals(compDeriv[i + 1], refDeriv[0][i], 1.1e-9);
}
}
......
......@@ -265,7 +265,7 @@ public class FieldViennaOneModelTest {
}
for (int i = 0; i < 6; i++) {
Assert.assertEquals(compDelay[i + 1], refDeriv[0][i], 1.2e-11);
Assert.assertEquals(compDelay[i + 1], refDeriv[0][i], 3.0e-11);
}
}
......
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment