-c altitude now works for FileEc

src/Calibrator/Altitude.cpp

@@ -30,6 +30,6 @@ bool CalibratorAltitude::calibrateCore(File& iFile, const ParameterFile* iParame
 
 std::string CalibratorAltitude::description() {
    std::stringstream ss;
-   ss << Util::formatDescription("-c altitude", "Changes the altitudes in the file to the altitudes in the parameter file.") << std::endl;
+   ss << Util::formatDescription("-c altitude", "Changes the altitudes in the file to the altitudes in the parameter file. If the file does not have an altitude field, then a new one is not created.") << std::endl;
    return ss.str();
 }

src/File/Ec.cpp

@@ -151,6 +151,7 @@ void FileEc::writeCore(std::vector<Variable::Type> iVariables) {
 
    writeTimes();
    writeReferenceTime();
+   writeAltitude();
    for(int v = 0; v < iVariables.size(); v++) {
       Variable::Type varType = iVariables[v];
       std::string variable = getVariableName(varType);
@@ -390,3 +391,57 @@ std::string FileEc::description() {
    ss << Util::formatDescription("type=ec", "ECMWF ensemble file") << std::endl;
    return ss.str();
 }
+
+void FileEc::writeAltitude() const {
+   if(!hasVar("altitude")) {
+      return;
+   }
+   int vElev = getVar("altitude");
+   int numDims = getNumDims(vElev);
+   if(numDims == 1) {
+      Util::error("Cannot write altitude when the variable only has one dimension");
+   }
+
+   int N = getNumDims(vElev);
+   long count[N];
+   int size = 1;
+   int indexLat = Util::MV;
+   int indexLon = Util::MV;
+   int dims[N];
+   nc_inq_vardimid(mFile, vElev, dims);
+   for(int i = 0; i < N; i++) {
+      if(dims[i] == getLatDim()) {
+         count[i] = getNumLat();
+         size *= count[i];
+         indexLat = i;
+      }
+      else if(dims[i] == getLonDim()) {
+         count[i] = getNumLon();
+         size *= count[i];
+         indexLon = i;
+      }
+      else {
+         count[i] = 1;
+      }
+   }
+   float MV = getMissingValue(vElev);
+   float* values = new float[size];
+   vec2 elevs = getElevs();
+   for(int i = 0; i < getNumLat(); i++) {
+      for(int j = 0; j < getNumLon(); j++) {
+         float elev = elevs[i][j];
+         if(Util::isValid(MV) && !Util::isValid(elev))
+            elev = MV;
+         // Latitude dimension is ordered first
+         if(indexLat < indexLon) {
+            values[i*getNumLon() + j] = elev;
+         }
+         // Longitude dimension is ordered first
+         else {
+            values[j*getNumLat() + i] = elev;
+         }
+      }
+   }
+   nc_put_var_float(mFile, vElev, values);
+   delete[] values;
+}

src/File/Ec.h

@@ -29,6 +29,7 @@ class FileEc : public FileNetcdf {
 
       std::vector<int> mTimes;
       vec2 getGridValues(int iVariable) const;
+      void writeAltitude() const;
       int getLatDim() const;
       int getLonDim() const;
       int getLatVar() const;

src/Testing/CalibratorAltitude.cpp

@@ -16,7 +16,7 @@ namespace {
          virtual void TearDown() {
          }
    };
-   TEST_F(TestCalibratorAltitude, 10x10) {
+   TEST_F(TestCalibratorAltitude, arome) {
       FileArome from("testing/files/10x10.nc");
       ParameterFileNetcdf par(Options("file=testing/files/10x10_param_zero_altitude.nc"));
       CalibratorAltitude cal = CalibratorAltitude(Options());
@@ -36,6 +36,26 @@ namespace {
       EXPECT_FLOAT_EQ(304, (*after)(5,9,0));
       EXPECT_FLOAT_EQ(320, (*after)(0,9,0));
    }
+   TEST_F(TestCalibratorAltitude, ec) {
+      FileEc from("testing/files/10x10_ec.nc");
+      ParameterFileNetcdf par(Options("file=testing/files/10x10_param_zero_altitude.nc"));
+      CalibratorAltitude cal = CalibratorAltitude(Options());
+
+      cal.calibrate(from, &par);
+
+      // Elevations should all be 0
+      vec2 elevs = from.getElevs();
+      for(int i = 0; i < from.getNumLat(); i++) {
+         for(int j = 0; j < from.getNumLon(); j++) {
+            EXPECT_FLOAT_EQ(0, elevs[i][j]);
+         }
+      }
+      // Shouldn't have changed anything else
+      FieldPtr after = from.getField(Variable::T, 0);
+      EXPECT_FLOAT_EQ(301, (*after)(5,2,0));
+      EXPECT_FLOAT_EQ(304, (*after)(5,9,0));
+      EXPECT_FLOAT_EQ(320, (*after)(0,9,0));
+   }
    // Should not calibrate when a parameter file with no locaions is used
    TEST_F(TestCalibratorAltitude, locationIndependent) {
       ::testing::FLAGS_gtest_death_test_style = "threadsafe";