Skip to content

Commit

Permalink
Cleanup of the of the ioda converter (#634)
Browse files Browse the repository at this point in the history
* wip

* code tidy, simplified nc reader

* added empty example

* addressed rmclaren s review

* override to final
  • Loading branch information
guillaumevernieres authored Sep 20, 2023
1 parent a028e29 commit 945c359
Show file tree
Hide file tree
Showing 9 changed files with 162 additions and 115 deletions.
108 changes: 84 additions & 24 deletions utils/obsproc/NetCDFToIodaConverter.h
Original file line number Diff line number Diff line change
Expand Up @@ -14,46 +14,106 @@
#include "oops/util/DateTime.h"
#include "oops/util/Duration.h"
#include "oops/util/Logger.h"
#include "oops/util/missingValues.h"

namespace gdasapp {
// A simple data structure to organize the info to provide to the ioda
// writter
struct IodaVars {
int location;
int nVars;
Eigen::ArrayXf obsVal;
Eigen::ArrayXf obsError;
Eigen::ArrayXi preQc;
std::string units; // reference date for epoch time
};

// Base class for the converters
class NetCDFToIodaConverter {
public:
// Constructor that stores the configuration as a data member
explicit NetCDFToIodaConverter(const eckit::Configuration & fullConfig) {
// time window info
std::string winbegin;
std::string winend;
fullConfig.get("window begin", winbegin);
fullConfig.get("window end", winend);
this->windowBegin_ = util::DateTime(winbegin);
this->windowEnd_ = util::DateTime(winend);

// get input netcdf files
fullConfig.get("input files", this->inputFilenames_);

// ioda output file name
this->outputFilename_ = fullConfig.getString("output file");
}

// Method to write out a IODA file
void WriteToIoda() {
// Constructor: Stores the configuration as a data members
explicit NetCDFToIodaConverter(const eckit::Configuration & fullConfig) {
// time window info
std::string winbegin;
std::string winend;
std::string obsvar;
fullConfig.get("window begin", winbegin);
fullConfig.get("window end", winend);
fullConfig.get("variable", obsvar);
windowBegin_ = util::DateTime(winbegin);
windowEnd_ = util::DateTime(winend);
variable_ = obsvar;
oops::Log::info() << "--- Window begin: " << winbegin << std::endl;
oops::Log::info() << "--- Window end: " << winend << std::endl;
oops::Log::info() << "--- Variable: " << obsvar << std::endl;

// get input netcdf files
fullConfig.get("input files", inputFilenames_);
oops::Log::info() << "--- Input files: " << inputFilenames_ << std::endl;

// ioda output file name
outputFilename_ = fullConfig.getString("output file");
oops::Log::info() << "--- Output files: " << outputFilename_ << std::endl;
}

// Method to write out a IODA file (writter called in destructor)
void writeToIoda() {
// Create empty group backed by HDF file
ioda::Group group =
ioda::Engines::HH::createFile(
this->outputFilename_,
outputFilename_,
ioda::Engines::BackendCreateModes::Truncate_If_Exists);
this->ReadFromNetCDF(group);

// Extract ioda variables from the provider's files
std::string fileName = inputFilenames_[0]; // TODO(Guillaume): make it work for a list
gdasapp::IodaVars iodaVars;
providerToIodaVars(fileName, iodaVars);
oops::Log::debug() << "--- iodaVars.location: " << iodaVars.location << std::endl;
oops::Log::debug() << "--- iodaVars.obsVal: " << iodaVars.obsVal << std::endl;

// Update the group with the location dimension
ioda::NewDimensionScales_t
newDims {ioda::NewDimensionScale<int>("Location", iodaVars.location)};
ioda::ObsGroup ogrp = ioda::ObsGroup::generate(group, newDims);

// Set up the float creation parameters
ioda::VariableCreationParameters float_params;
float_params.chunk = true; // allow chunking
float_params.compressWithGZIP(); // compress using gzip
float missing_value = util::missingValue(missing_value);
float_params.setFillValue<float>(missing_value);

// Create the IODA variables
ioda::Variable adtIodaDatetime =
ogrp.vars.createWithScales<float>("Metadata/dateTime",
{ogrp.vars["Location"]}, float_params);
// TODO(Mindo): Get the date info from the netcdf file
adtIodaDatetime.atts.add<std::string>("units", {"seconds since 9999-04-15T12:00:00Z"}, {1});

ioda::Variable adtIodaObsVal =
ogrp.vars.createWithScales<float>("ObsValue/"+variable_,
{ogrp.vars["Location"]}, float_params);
ioda::Variable adtIodaObsErr =
ogrp.vars.createWithScales<float>("ObsError/"+variable_,
{ogrp.vars["Location"]}, float_params);

// Write adt obs value to group
adtIodaObsVal.writeWithEigenRegular(iodaVars.obsVal);

// Write adt obs error to group
adtIodaObsErr.writeWithEigenRegular(iodaVars.obsError);
}

private:
// Virtual method to read a netcdf file and store the relevant info in a group
virtual void ReadFromNetCDF(ioda::Group &) = 0;

// Virtual method that reads the provider's netcdf file and store the relevant
// info in a IodaVars struct
virtual void providerToIodaVars(const std::string fileName, gdasapp::IodaVars & iodaVars) = 0;

protected:
util::DateTime windowBegin_;
util::DateTime windowEnd_;
std::vector<std::string> inputFilenames_;
std::string outputFilename_;
std::string variable_;
};
} // namespace gdasapp
68 changes: 21 additions & 47 deletions utils/obsproc/Rads2Ioda.h
Original file line number Diff line number Diff line change
Expand Up @@ -13,75 +13,49 @@

#include "NetCDFToIodaConverter.h"

#include "oops/util/missingValues.h"

namespace gdasapp {

class Rads2Ioda : public NetCDFToIodaConverter {
public:
explicit Rads2Ioda(const eckit::Configuration & fullConfig)
: NetCDFToIodaConverter(fullConfig) {
oops::Log::info() << "Window begin: " << this->windowBegin_ << std::endl;
oops::Log::info() << "Window end: " << this->windowEnd_ << std::endl;
oops::Log::info() << "IODA output file name: " << this->outputFilename_ << std::endl;
}

// Read netcdf file and populate group
void ReadFromNetCDF(ioda::Group & group) override {
std::string fileName = this->inputFilenames_[0]; // TODO(Guillaume): make it work for a list
// Read netcdf file and populate iodaVars
void providerToIodaVars(const std::string fileName, gdasapp::IodaVars & iodaVars) final {
oops::Log::info() << "Processing files provided by the RADS" << std::endl;

// Set nVars to 1
iodaVars.nVars = 1;

// Open the NetCDF file in read-only mode
netCDF::NcFile ncFile(fileName, netCDF::NcFile::read);

// Get dimensions
int location = ncFile.getDim("time").getSize();
int nVars = 1;
iodaVars.location = ncFile.getDim("time").getSize();
oops::Log::debug() << "--- iodaVars.location: " << iodaVars.location << std::endl;

// Allocate memory
iodaVars.obsVal = Eigen::ArrayXf(iodaVars.location);
iodaVars.obsError = Eigen::ArrayXf(iodaVars.location);
iodaVars.preQc = Eigen::ArrayXi(iodaVars.location);

// Get adt_egm2008 obs values and attributes
netCDF::NcVar adtNcVar = ncFile.getVar("adt_egm2008");
int adtObsVal[location]; // NOLINT (can't pass vector to getVar below)
int adtObsVal[iodaVars.location]; // NOLINT (can't pass vector to getVar below)
adtNcVar.getVar(adtObsVal);
std::string units;
adtNcVar.getAtt("units").getValues(units);
float scaleFactor;
adtNcVar.getAtt("scale_factor").getValues(&scaleFactor);

// Update the group with 1 dimension: location
ioda::NewDimensionScales_t newDims {ioda::NewDimensionScale<int>("Location", location)};
ioda::ObsGroup ogrp = ioda::ObsGroup::generate(group, newDims);

// Set up the float creation parameters
ioda::VariableCreationParameters float_params;
float_params.chunk = true; // allow chunking
float_params.compressWithGZIP(); // compress using gzip
float missing_value = util::missingValue(missing_value);
float_params.setFillValue<float>(missing_value);

// Create the IODA variables
ioda::Variable adtIodaDatetime =
ogrp.vars.createWithScales<float>("Metadata/dateTime",
{ogrp.vars["Location"]}, float_params);
// TODO(Mindo): Get the date info from the netcdf file
adtIodaDatetime.atts.add<std::string>("units", {"seconds since 9999-04-15T12:00:00Z"}, {1});

ioda::Variable adtIodaObsVal =
ogrp.vars.createWithScales<float>("ObsValue/absoluteDynamicTopography",
{ogrp.vars["Location"]}, float_params);
ioda::Variable adtIodaObsErr =
ogrp.vars.createWithScales<float>("ObsError/absoluteDynamicTopography",
{ogrp.vars["Location"]}, float_params);

// Write adt obs value to group
Eigen::ArrayXf tmpvar(location);
for (int i = 0; i <= location; i++) {
tmpvar[i] = static_cast<float>(adtObsVal[i])*scaleFactor;
for (int i = 0; i <= iodaVars.location; i++) {
iodaVars.obsVal(i) = static_cast<float>(adtObsVal[i])*scaleFactor;
}
adtIodaObsVal.writeWithEigenRegular(tmpvar);

// Write adt obs error to group
for (int i = 0; i <= location; i++) {
tmpvar[i] = 0.1;
// Do something for obs error
for (int i = 0; i <= iodaVars.location; i++) {
iodaVars.obsError(i) = 0.1;
}
adtIodaObsErr.writeWithEigenRegular(tmpvar);
};
};
}; // class Rads2Ioda
} // namespace gdasapp
10 changes: 5 additions & 5 deletions utils/obsproc/applications/CMakeLists.txt
Original file line number Diff line number Diff line change
@@ -1,9 +1,9 @@
list( APPEND gdasapp_rads2ioda_src_files
gdas_rads2ioda.cc
gdas_rads2ioda.h
gdas_obsprovider2ioda.cc
gdas_obsprovider2ioda.h
)
ecbuild_add_executable( TARGET gdas_rads2ioda.x
ecbuild_add_executable( TARGET gdas_obsprovider2ioda.x
SOURCES ${gdasapp_rads2ioda_src_files} )

target_compile_features( gdas_rads2ioda.x PUBLIC cxx_std_17)
target_link_libraries( gdas_rads2ioda.x PUBLIC oops ioda NetCDF::NetCDF_CXX)
target_compile_features( gdas_obsprovider2ioda.x PUBLIC cxx_std_17)
target_link_libraries( gdas_obsprovider2ioda.x PUBLIC oops ioda NetCDF::NetCDF_CXX)
8 changes: 8 additions & 0 deletions utils/obsproc/applications/gdas_obsprovider2ioda.cc
Original file line number Diff line number Diff line change
@@ -0,0 +1,8 @@
#include "gdas_obsprovider2ioda.h"
#include "oops/runs/Run.h"

int main(int argc, char ** argv) {
oops::Run run(argc, argv);
gdasapp::ObsProvider2IodaApp obsprovider2ioda;
return run.execute(obsprovider2ioda);
}
41 changes: 41 additions & 0 deletions utils/obsproc/applications/gdas_obsprovider2ioda.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,41 @@
#pragma once

#include <string>

#include "eckit/config/LocalConfiguration.h"
#include "oops/mpi/mpi.h"
#include "oops/runs/Application.h"

#include "../Rads2Ioda.h"

namespace gdasapp {
class ObsProvider2IodaApp : public oops::Application {
public:
explicit ObsProvider2IodaApp(const eckit::mpi::Comm & comm = oops::mpi::world())
: Application(comm) {}
static const std::string classname() {return "gdasapp::ObsProvider2IodaApp";}

int execute(const eckit::Configuration & fullConfig, bool /*validate*/) const {
// Get the file provider string identifier from the config
std::string provider;
fullConfig.get("provider", provider);

if (provider == "RADS") {
Rads2Ioda conv2ioda(fullConfig);
conv2ioda.writeToIoda();
} else if (provider == "GHRSST") {
oops::Log::info() << "Comming soon!" << std::endl;
} else {
oops::Log::info() << "Provider not implemented" << std::endl;
return 1;
}
return 0;
}
// -----------------------------------------------------------------------------
private:
std::string appname() const {
return "gdasapp::ObsProvider2IodaApp";
}
// -----------------------------------------------------------------------------
};
} // namespace gdasapp
8 changes: 0 additions & 8 deletions utils/obsproc/applications/gdas_rads2ioda.cc

This file was deleted.

30 changes: 0 additions & 30 deletions utils/obsproc/applications/gdas_rads2ioda.h

This file was deleted.

2 changes: 1 addition & 1 deletion utils/test/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@ ecbuild_add_test( TARGET test_gdasapp_util_prepdata

# Test the RADS to IODA converter
ecbuild_add_test( TARGET test_gdasapp_util_rads2ioda
COMMAND ${CMAKE_BINARY_DIR}/bin/gdas_rads2ioda.x
COMMAND ${CMAKE_BINARY_DIR}/bin/gdas_obsprovider2ioda.x
ARGS "../testinput/gdas_rads2ioda.yaml"
LIBS gdas-utils
WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR}/obsproc)
2 changes: 2 additions & 0 deletions utils/test/testinput/gdas_rads2ioda.yaml
Original file line number Diff line number Diff line change
@@ -1,5 +1,7 @@
provider: RADS
window begin: 2018-04-15T06:00:00Z
window end: 2018-04-15T12:00:00Z
variable: absoluteDynamicTopography
output file: rads_adt_j3_2021182.ioda.nc
input files:
- rads_adt_j3_2021182.nc4

0 comments on commit 945c359

Please sign in to comment.