WSL/SLF GitLab Repository

Commit 3555d63d authored by Mathias Bavay's avatar Mathias Bavay
Browse files

Small optimization into Sun: if setting the lat/lon to the same values as...

Small optimization into Sun: if setting the lat/lon to the same values as internally stored, nothing is changed (so it does not costs anything). 

The CNRMIO pluggin has been changed by removing the grids support (this should be performed by the generic NetCDF plugin) and adding what is necessary for Crocus to be able to use the generated meteo time series.
parent 6a51065f
......@@ -75,6 +75,14 @@ void SunObject::setLatLon(const double& i_latitude, const double& i_longitude, c
throw NoDataException(ss.str(), AT);
}
if (i_latitude==latitude && i_longitude==longitude) {
if (altitude==i_altitude) return; //everything is the same, nothing to do
altitude = i_altitude;
beam_toa = beam_direct = beam_diffuse = IOUtils::nodata;
return;
}
position.reset();
latitude = i_latitude;
longitude = i_longitude;
......
......@@ -20,6 +20,7 @@
#include <meteoio/meteoStats/libinterpol1D.h>
#include <meteoio/Timer.h>
#include <meteoio/MathOptim.h>
#include <meteoio/meteoLaws/Atmosphere.h>
#include <meteoio/plugins/libncpp.h>
#include <cmath>
......@@ -70,19 +71,12 @@ namespace mio {
*/
const double CNRMIO::plugin_nodata = -9999999.; //CNRM-GAME nodata value
const double CNRMIO::epsilon = 1.0e-10; //when comparing timestamps
const std::string CNRMIO::cf_time = "time";
const std::string CNRMIO::cf_units = "units";
const std::string CNRMIO::cf_days = "days since ";
const std::string CNRMIO::cf_hours = "hours since ";
const std::string CNRMIO::cf_seconds = "seconds since ";
const std::string CNRMIO::cf_latitude = "lat";
const std::string CNRMIO::cf_longitude = "lon";
const std::string CNRMIO::cf_altitude = "z";
const std::string CNRMIO::cf_ta = "temperature";
const std::string CNRMIO::cf_rh = "humidity";
const std::string CNRMIO::cf_p = "pressure";
const std::string CNRMIO::cnrm_points = "Number_of_points";
const std::string CNRMIO::cnrm_latitude = "LAT";
......@@ -90,17 +84,17 @@ const std::string CNRMIO::cnrm_longitude = "LON";
const std::string CNRMIO::cnrm_altitude = "ZS";
const std::string CNRMIO::cnrm_aspect = "aspect";
const std::string CNRMIO::cnrm_slope = "slope";
const std::string CNRMIO::cnrm_uref = "uref";
const std::string CNRMIO::cnrm_zref = "zref";
const std::string CNRMIO::cnrm_uref = "UREF";
const std::string CNRMIO::cnrm_zref = "ZREF";
const std::string CNRMIO::cnrm_ta = "Tair";
const std::string CNRMIO::cnrm_rh = "HUMREL";
const std::string CNRMIO::cnrm_vw = "Wind";
const std::string CNRMIO::cnrm_dw = "Wind_DIR";
const std::string CNRMIO::cnrm_qair = "Qair";
const std::string CNRMIO::cnrm_co2air = "CO2air";
const std::string CNRMIO::cnrm_theorsw = "theorSW";
const std::string CNRMIO::cnrm_iswr = "theorSW";
const std::string CNRMIO::cnrm_neb = "NEB";
const std::string CNRMIO::cnrm_psum = "Rainf";
const std::string CNRMIO::cnrm_rainf = "Rainf";
const std::string CNRMIO::cnrm_snowf = "Snowf";
const std::string CNRMIO::cnrm_swr_direct = "DIR_SWdown";
const std::string CNRMIO::cnrm_swr_diffuse = "SCA_SWdown";
......@@ -119,11 +113,11 @@ bool CNRMIO::initStaticData()
paramname[cnrm_qair] = IOUtils::npos; // not a standard MeteoIO parameter
paramname[cnrm_co2air] = IOUtils::npos; // not a standard MeteoIO parameter
paramname[cnrm_neb] = IOUtils::npos; // not a standard MeteoIO parameter
paramname[cnrm_theorsw] = IOUtils::npos; // not a standard MeteoIO parameter
paramname[cnrm_iswr] = MeteoData::ISWR;
paramname[cnrm_rh] = MeteoData::RH;
paramname[cnrm_vw] = MeteoData::VW;
paramname[cnrm_dw] = MeteoData::DW;
paramname[cnrm_psum] = IOUtils::npos;
paramname[cnrm_rainf] = IOUtils::npos;
paramname[cnrm_snowf] = IOUtils::npos;
paramname[cnrm_swr_direct] = IOUtils::npos;
paramname[cnrm_swr_diffuse] = IOUtils::npos;
......@@ -136,11 +130,10 @@ bool CNRMIO::initStaticData()
map_name["P"] = cnrm_p;
map_name["VW"] = cnrm_vw;
map_name["DW"] = cnrm_dw;
map_name["ISWR"] = cnrm_swr_direct;
map_name["PSUM"] = cnrm_psum;
map_name["ISWR"] = cnrm_iswr;
map_name["PSUM"] = cnrm_rainf;
map_name[cnrm_co2air] = cnrm_co2air;
map_name[cnrm_qair] = cnrm_qair;
map_name[cnrm_theorsw] = cnrm_theorsw;
map_name[cnrm_neb] = cnrm_neb;
return true;
......@@ -174,102 +167,6 @@ void CNRMIO::parseInputOutputSection()
cfg.getValue("ZREF", "Output", zref, IOUtils::nothrow);
}
void CNRMIO::read2DGrid(Grid2DObject& grid_out, const std::string& arguments)
{
vector<string> vec_argument;
IOUtils::readLineToVec(arguments, vec_argument, ':');
if (vec_argument.size() == 2) {
read2DGrid_internal(grid_out, vec_argument[0], vec_argument[1]);
} else {
throw InvalidArgumentException("The format for the arguments to CNRMIO::read2DGrid is filename:varname", AT);
}
}
void CNRMIO::read2DGrid(Grid2DObject& grid_out, const MeteoGrids::Parameters& parameter, const Date& date)
{
const string filename = cfg.get("GRID2DFILE", "Input");
const string varname( get_varname(parameter) );
read2DGrid_internal(grid_out, filename, varname, date);
}
void CNRMIO::read2DGrid_internal(Grid2DObject& grid_out, const std::string& filename, const std::string& varname, const Date& date)
{
const bool is_record = (date != Date());
size_t lat_index = 0, lon_index = 1;
if (!FileUtils::fileExists(filename)) throw AccessException(filename, AT); //prevent invalid filenames
int ncid, varid;
vector<int> dimid, dim_varid;
ncpp::open_file(filename, NC_NOWRITE, ncid);
ncpp::get_variable(ncid, varname, varid);
vector<string> dimname;
vector<size_t> dimlen;
ncpp::get_dimension(ncid, varname, varid, dimid, dim_varid, dimname, dimlen);
if (is_record) { // In case we're reading a record the first index is always the record index
lat_index = 1;
lon_index = 2;
if (dimid.size()!=3 || dimlen[0]<1 || dimlen[lat_index]<2 || dimlen[lon_index]<2)
throw IOException("Variable '" + varname + "' may only have three dimensions, all have to at least have length 1", AT);
} else if (dimid.size()==3 && dimlen[0]==1) { //in case the variable is associated with a 1 element time dimension
lat_index = 1;
lon_index = 2;
if (dimlen[lat_index]<2 || dimlen[lon_index]<2)
throw IOException("All dimensions for variable '" + varname + "' have to at least have length 1", AT);
} else if (dimid.size()!=2 || dimlen[lat_index]<2 || dimlen[lon_index]<2) {
throw IOException("Variable '" + varname + "' may only have two dimensions and both have to have length >1", AT);
}
double *lat = new double[dimlen[lat_index]];
double *lon = new double[dimlen[lon_index]];
double *grid = new double[dimlen[lat_index]*dimlen[lon_index]];
ncpp::read_data(ncid, dimname[lat_index], dim_varid[lat_index], lat);
ncpp::read_data(ncid, dimname[lon_index], dim_varid[lon_index], lon);
if (is_record) {
const size_t pos = ncpp::find_record(ncid, CNRMIO::cf_time, dimid[0], date.getModifiedJulianDate());
if (pos == IOUtils::npos)
throw IOException("No record for date " + date.toString(Date::ISO), AT);
ncpp::read_data(ncid, varname, varid, pos, dimlen[lat_index], dimlen[lon_index], grid);
} else {
ncpp::read_data(ncid, varname, varid, grid);
}
double missing_value = plugin_nodata;
if (ncpp::check_attribute(ncid, varid, "missing_value"))
ncpp::get_attribute(ncid, varname, varid, "missing_value", missing_value);
ncpp::copy_grid(coordin, coordinparam, dimlen[lat_index], dimlen[lon_index], lat, lon, grid, missing_value, grid_out);
//handle data packing if necessary
if (ncpp::check_attribute(ncid, varid, "scale_factor")) {
double scale_factor = 1.;
ncpp::get_attribute(ncid, varname, varid, "scale_factor", scale_factor);
grid_out.grid2D *= scale_factor;
}
if (ncpp::check_attribute(ncid, varid, "add_offset")) {
double add_offset = 0.;
ncpp::get_attribute(ncid, varname, varid, "add_offset", add_offset);
grid_out.grid2D += add_offset;
}
ncpp::close_file(filename, ncid);
delete[] lat; delete[] lon; delete[] grid;
}
void CNRMIO::readDEM(DEMObject& dem_out)
{
const string filename = cfg.get("DEMFILE", "Input");
const string varname = cfg.get("DEMVAR", "Input");
read2DGrid_internal(dem_out, filename, varname);
}
void CNRMIO::readStationData(const Date&, std::vector<StationData>& vecStation)
{
if (!vecMetaData.empty()) { // We already have meta data
......@@ -419,7 +316,7 @@ void CNRMIO::copy_data(const int& ncid, const std::map<std::string, size_t>& map
const size_t param = map_parameters.find(varname)->second; //must exist, at this point we know it does
if (param == IOUtils::npos) {
if ((varname == cnrm_snowf) || (varname == cnrm_psum)) {
if ((varname == cnrm_snowf) || (varname == cnrm_rainf)) {
int varid;
ncpp::get_variable(ncid, cnrm_timestep, varid);
ncpp::read_value(ncid, cnrm_timestep, varid, multiplier);
......@@ -503,7 +400,7 @@ void CNRMIO::get_parameters(const int& ncid, std::map<std::string, size_t>& map_
if (name==cnrm_qair)
index = meteo_data.addParameter( "SH" );
if ((name == cnrm_theorsw) || (name == cnrm_co2air) || (name == cnrm_neb)) {
if ((name == cnrm_swr_diffuse) || (name == cnrm_swr_direct) || (name == cnrm_co2air) || (name == cnrm_neb)) {
index = meteo_data.addParameter( name );
}
......@@ -637,7 +534,6 @@ void CNRMIO::writeMeteoData(const std::vector< std::vector<MeteoData> >& vecMete
const string file_and_path = path + "/" + filename;
int ncid, did_time, vid_time, did_points;
bool create_time = false, create_points = false, create_locations = false, create_variables = false;
const bool exists = FileUtils::fileExists(file_and_path);
if (exists) remove(file_and_path.c_str()); // NOTE: file is deleted if it exists
......@@ -651,7 +547,7 @@ void CNRMIO::writeMeteoData(const std::vector< std::vector<MeteoData> >& vecMete
map_data[cnrm_uref] = new double[number_of_stations];
map_data[cnrm_zref] = new double[number_of_stations];
//compute data set start julian date
//compute data set start julian date and add Qair
Date ref_date;
for (size_t ii=0; ii<number_of_stations; ii++) {
if (vecMeteo[ii].empty()) continue;
......@@ -666,13 +562,10 @@ void CNRMIO::writeMeteoData(const std::vector< std::vector<MeteoData> >& vecMete
get_parameters(ref_date.getJulian(), vecMeteo, map_param_name, map_data, dates);
ncpp::create_file(file_and_path, NC_CLASSIC_MODEL, ncid);
create_time = create_points = create_locations = create_variables = true;
if (create_time) create_time_dimension(ref_date, ncid, did_time, vid_time);
if (create_points) ncpp::add_dimension(ncid, cnrm_points, number_of_stations, did_points);
if (create_locations) create_meta_data(ncid, did_points, map_data, varid);
if (create_variables) create_parameters(ncid, did_time, did_points, number_of_records, number_of_stations, map_param_name, map_data, varid);
create_time_dimension(ref_date, ncid, did_time, vid_time);
ncpp::add_dimension(ncid, cnrm_points, number_of_stations, did_points);
create_meta_data(ncid, did_points, map_data, varid);
create_parameters(ncid, did_time, did_points, number_of_records, number_of_stations, map_param_name, map_data, varid);
ncpp::end_definitions(file_and_path, ncid);
copy_data(number_of_stations, number_of_records, vecMeteo, map_param_name, map_data);
......@@ -698,32 +591,70 @@ void CNRMIO::copy_data(const size_t& number_of_stations, const size_t& number_of
for (map<size_t, string>::const_iterator it = map_param_name.begin(); it != map_param_name.end(); ++it) {
const size_t param( it->first );
const string varname( it->second );
bool simple_copy = false, multiply_copy = false;
double multiplier = IOUtils::nodata;
double* data = map_data_2D[varname];
if (param == MeteoData::RH) {
multiplier = 100.;
multiply_copy = true;
} else if (param == MeteoData::PSUM) {
multiply_copy = true;
multiplier = 1./3600.;
} else {
simple_copy = true;
if (varname==cnrm_qair) { //special processing for Qair
for (size_t ii=0; ii<number_of_stations; ++ii) {
const double altitude = vecMeteo[ii].front().meta.position.getAltitude();
for (size_t jj=0; jj<number_of_records; ++jj) {
const double TA = vecMeteo[ii][jj](MeteoData::TA);
const double RH = vecMeteo[ii][jj](MeteoData::RH);
if (altitude!=IOUtils::nodata && TA!=IOUtils::nodata && RH!=IOUtils::nodata)
data[jj*number_of_stations + ii] = Atmosphere::relToSpecHumidity(altitude, TA, RH);
else
data[jj*number_of_stations + ii] = plugin_nodata;
}
}
continue;
}
if (varname==cnrm_rainf) {
for (size_t ii=0; ii<number_of_stations; ++ii) {
if (number_of_records>0) {
map_data_2D[cnrm_rainf][ii] = plugin_nodata;
map_data_2D[cnrm_snowf][ii] = plugin_nodata;
}
for (size_t jj=1; jj<number_of_records; ++jj) {
const double acc_period = (vecMeteo[ii][jj].date.getJulian() - vecMeteo[ii][jj-1].date.getJulian()) * (24.*3600.); //in seconds
if (acc_period==0.) continue; //this should not happen, but...
const double psum = vecMeteo[ii][jj](MeteoData::PSUM);
const double phase = vecMeteo[ii][jj](MeteoData::PSUM_PH);
if (psum!=IOUtils::nodata) {
if (phase!=IOUtils::nodata) {
map_data_2D[cnrm_rainf][jj*number_of_stations + ii] = psum * phase / acc_period;
map_data_2D[cnrm_snowf][jj*number_of_stations + ii] = psum * (1.-phase) / acc_period;
} else {
map_data_2D[cnrm_rainf][jj*number_of_stations + ii] = psum / acc_period;
map_data_2D[cnrm_snowf][jj*number_of_stations + ii] = plugin_nodata;
}
} else {
map_data_2D[cnrm_rainf][jj*number_of_stations + ii] = plugin_nodata;
map_data_2D[cnrm_snowf][jj*number_of_stations + ii] = plugin_nodata;
}
}
}
continue;
}
if (varname==cnrm_snowf) continue; //this is handled by cnrm_rainf
for (size_t ii=0; ii<number_of_stations; ++ii) {
const size_t nr_of_parameters = vecMeteo[ii].front().getNrOfParameters();
if (param>=nr_of_parameters) { //unknown parameters are filled with nodata
for (size_t jj=0; jj<number_of_records; ++jj)
data[jj*number_of_stations + ii] = plugin_nodata;
continue;
}
for (size_t jj=0; jj<number_of_records; ++jj) {
const double& value = vecMeteo[ii][jj](param);
if (value == IOUtils::nodata) {
data[jj*number_of_stations + ii] = plugin_nodata;
} else if (simple_copy) {
} else {
data[jj*number_of_stations + ii] = value;
} else if (multiply_copy) {
data[jj*number_of_stations + ii] = value * multiplier;
}
}
}
......@@ -810,29 +741,34 @@ void CNRMIO::get_parameters(const double& ref_julian, const std::vector< std::ve
}
const size_t nr_of_parameters = (!vecMeteo[0].empty())? vecMeteo[0][0].getNrOfParameters() : 0 ;
vector<bool> vec_param_in_use(nr_of_parameters, false);
vector<string> vec_param_name(nr_of_parameters, "");
//Check consistency, dates must be existent everywhere
bool inconsistent = false;
for (size_t ii=0; ii<vecMeteo.size(); ++ii) {
if (vecMeteo[ii].empty()) continue;
//add metadata
const StationData sd( vecMeteo[ii].front().meta );
const double slope = sd.getSlopeAngle();
const double azi = sd.getAzimuth();
map_data_1D[cnrm_latitude][ii] = toNetcdfNodata( sd.position.getLat() );
map_data_1D[cnrm_longitude][ii] = toNetcdfNodata( sd.position.getLon() );
map_data_1D[cnrm_altitude][ii] = toNetcdfNodata( sd.position.getAltitude() );
map_data_1D[cnrm_slope][ii] = (slope==IOUtils::nodata)? 0 : slope;
map_data_1D[cnrm_aspect][ii] = (azi==IOUtils::nodata)? 0 : azi;
map_data_1D[cnrm_uref][ii] = uref;
map_data_1D[cnrm_zref][ii] = zref;
if (number_of_records != vecMeteo[ii].size()) inconsistent = true;
for (size_t jj=0; jj<vecMeteo[ii].size(); ++jj) {
const MeteoData& meteo_data = vecMeteo[ii][jj];
const MeteoData& meteo_data( vecMeteo[ii][jj] );
//check time steps consistency compared to what we declared in the metadata
const int local_date = static_cast<int>( Optim::round( (meteo_data.date.getJulian() - ref_julian) * 24.*3600. ) );
if (!IOUtils::checkEpsilonEquality(dates[jj], local_date, CNRMIO::epsilon))
inconsistent = true;
if (jj == 0) {
map_data_1D[cnrm_latitude][ii] = toNetcdfNodata( meteo_data.meta.position.getLat() );
map_data_1D[cnrm_longitude][ii] = toNetcdfNodata( meteo_data.meta.position.getLon() );
map_data_1D[cnrm_altitude][ii] = toNetcdfNodata( meteo_data.meta.position.getAltitude() );
map_data_1D[cnrm_slope][ii] = toNetcdfNodata( meteo_data.meta.getSlopeAngle() );
map_data_1D[cnrm_aspect][ii] = toNetcdfNodata( meteo_data.meta.getAzimuth() );
map_data_1D[cnrm_uref][ii] = uref;
map_data_1D[cnrm_zref][ii] = zref;
}
if (dates[jj] != local_date) inconsistent = true; //dates are ints in seconds
//Check which parameters are in use
for (size_t kk=0; kk<nr_of_parameters; ++kk) {
......@@ -852,6 +788,11 @@ void CNRMIO::get_parameters(const double& ref_julian, const std::vector< std::ve
if (vec_param_in_use[kk])
map_param_name[kk] = vec_param_name[kk];
}
map_param_name[nr_of_parameters] = cnrm_snowf; //so we can handle the precipitation splitting
map_param_name[nr_of_parameters+1] = cnrm_qair; //add Qair
map_param_name[nr_of_parameters+2] = cnrm_co2air; //this one needs to be there for Safran
map_param_name[nr_of_parameters+3] = cnrm_swr_diffuse; //this one needs to be there for Safran
map_param_name[nr_of_parameters+4] = cnrm_swr_direct; //this one needs to be there for Safran
double interval = (number_of_records>1)? static_cast<double>(dates[1] - dates[0]) : 0; //we consider this gives us the sampling rate
double* timestep = new double[1];
......@@ -859,132 +800,6 @@ void CNRMIO::get_parameters(const double& ref_julian, const std::vector< std::ve
map_data_1D[cnrm_timestep] = timestep;
}
void CNRMIO::write2DGrid(const Grid2DObject& grid_in, const std::string& arguments)
{
// arguments is a string of the format filname:varname
vector<string> vec_argument;
IOUtils::readLineToVec(arguments, vec_argument, ':');
if (vec_argument.size() != 2)
throw InvalidArgumentException("The format for the arguments to CNRMIO::write2DGrid is filename:varname", AT);
write2DGrid_internal(grid_in, vec_argument[0], vec_argument[1]);
}
void CNRMIO::write2DGrid(const Grid2DObject& grid_in, const MeteoGrids::Parameters& parameter, const Date& date)
{
const string filename = cfg.get("GRID2DFILE", "Output");
const string varname = get_varname(parameter);
write2DGrid_internal(grid_in, filename, varname, date);
}
void CNRMIO::write2DGrid_internal(const Grid2DObject& grid_in, const std::string& filename, const std::string& varname, const Date& date)
{
const bool is_record = (date != Date());
const bool exists = FileUtils::fileExists(filename);
double *lat_array = new double[grid_in.getNy()];
double *lon_array = new double[grid_in.getNx()];
double *data = new double[grid_in.getNy() * grid_in.getNx()];
ncpp::calculate_dimensions(grid_in, lat_array, lon_array);
ncpp::fill_grid_data(grid_in, data);
int ncid, did_lat, did_lon, did_time, vid_lat, vid_lon, vid_var, vid_time;
bool create_dimensions(false), create_variable(false), create_time(false);
if (exists) {
ncpp::open_file(filename, NC_WRITE, ncid);
//check of lat/lon are defined and consistent
if (ncpp::check_dim_var(ncid, cf_latitude) && ncpp::check_dim_var(ncid, cf_longitude)) {
check_consistency(ncid, grid_in, lat_array, lon_array, did_lat, did_lon, vid_lat, vid_lon);
} else {
create_dimensions = true;
}
if (is_record) {
//check if a time dimension/variable already exists
if (ncpp::check_dim_var(ncid, CNRMIO::cf_time)) {
ncpp::get_dimension(ncid, CNRMIO::cf_time, did_time);
ncpp::get_variable(ncid, CNRMIO::cf_time, vid_time);
} else {
create_time = true;
}
}
if (ncpp::check_variable(ncid, varname)) { // variable exists
ncpp::get_variable(ncid, varname, vid_var);
vector<int> dimid, dim_varid;
vector<string> dimname;
vector<size_t> dimlen;
ncpp::get_dimension(ncid, varname, vid_var, dimid, dim_varid, dimname, dimlen);
if (is_record) {
if ((dimname.size() != 3) || (dimname[0] != cf_time) || (dimname[1] != cf_latitude) || (dimname[2] != cf_longitude) || (dimlen[1]!=grid_in.getNy()) || (dimlen[2]!=grid_in.getNx()))
throw IOException("Variable '" + varname + "' already defined with different dimensions in file '"+ filename +"'", AT);
} else {
if ((dimname[0] != cf_latitude) || (dimname[1] != cf_longitude) || (dimlen[0]!=grid_in.getNy()) || (dimlen[1]!=grid_in.getNx()))
throw IOException("Variable '" + varname + "' already defined with different dimensions in file '"+ filename +"'", AT);
}
} else {
create_variable = true;
}
ncpp::start_definitions(filename, ncid);
} else {
if (!FileUtils::validFileAndPath(filename)) throw InvalidNameException(filename,AT);
ncpp::create_file(filename, NC_CLASSIC_MODEL, ncid);
ncpp::add_attribute(ncid, NC_GLOBAL, "Conventions", "CF-1.3");
create_variable = create_dimensions = true;
if (is_record) create_time = true;
}
if (create_dimensions) create_latlon_dimensions(ncid, grid_in, did_lat, did_lon, vid_lat, vid_lon);
if (create_time) create_time_dimension(date, ncid, did_time, vid_time);
if (is_record && create_variable) {
ncpp::add_3D_variable(ncid, varname, NC_DOUBLE, did_time, did_lat, did_lon, vid_var);
add_attributes_for_variable(ncid, vid_var, varname);
} else if (create_variable) {
ncpp::add_2D_variable(ncid, varname, NC_DOUBLE, did_lat, did_lon, vid_var);
add_attributes_for_variable(ncid, vid_var, varname);
}
ncpp::end_definitions(filename, ncid);
if (create_dimensions) {
ncpp::write_data(ncid, cf_latitude, vid_lat, lat_array);
ncpp::write_data(ncid, cf_longitude, vid_lon, lon_array);
}
if (is_record) {
const size_t pos_start = ncpp::add_record(ncid, CNRMIO::cf_time, vid_time, date.getModifiedJulianDate());
ncpp::write_data(ncid, varname, vid_var, grid_in.getNy(), grid_in.getNx(), pos_start, data);
} else {
ncpp::write_data(ncid, varname, vid_var, data);
}
ncpp::close_file(filename, ncid);
delete[] lat_array; delete[] lon_array; delete[] data;
}
void CNRMIO::create_latlon_dimensions(const int& ncid, const Grid2DObject& grid_in, int& did_lat, int& did_lon, int& vid_lat, int& vid_lon)
{
ncpp::add_dimension(ncid, cf_latitude, grid_in.getNy(), did_lat);
ncpp::add_1D_variable(ncid, cf_latitude, NC_DOUBLE, did_lat, vid_lat);
add_attributes_for_variable(ncid, vid_lat, cf_latitude);
ncpp::add_dimension(ncid, cf_longitude, grid_in.getNx(), did_lon);
ncpp::add_1D_variable(ncid, cf_longitude, NC_DOUBLE, did_lon, vid_lon);
add_attributes_for_variable(ncid, vid_lon, cf_longitude);
}
void CNRMIO::create_time_dimension(const Date& ref_date, const int& ncid, int& did_time, int& vid_time)
{
ncpp::add_dimension(ncid, CNRMIO::cf_time, NC_UNLIMITED, did_time);
......@@ -997,63 +812,29 @@ void CNRMIO::create_time_dimension(const Date& ref_date, const int& ncid, int& d
ncpp::add_attribute(ncid, vid_time, "units", "seconds since "+ref_string);
}
// When reading or writing gridded variables we should have a consistent naming
// scheme: http://cfconventions.org/1.6.html
std::string CNRMIO::get_varname(const MeteoGrids::Parameters& parameter)
{
string varname = MeteoGrids::getParameterName(parameter);
if (parameter == MeteoGrids::TA) varname = cnrm_ta;
else if (parameter == MeteoGrids::RH) varname = cnrm_rh;
else if (parameter == MeteoGrids::DEM) varname = cnrm_altitude;
else if (parameter == MeteoGrids::P) varname = cnrm_p;
else if (parameter == MeteoGrids::VW) varname = cnrm_vw;
else if (parameter == MeteoGrids::DW) varname = cnrm_dw;
else if (parameter == MeteoGrids::ILWR) varname = cnrm_ilwr;
else if (parameter == MeteoGrids::PSUM) varname = cnrm_psum; //HACK this should add snowf!
else if (parameter == MeteoGrids::SLOPE) varname = cnrm_slope;
else if (parameter == MeteoGrids::AZI) varname = cnrm_aspect;
//TODO: iswr=dir+diff
//TODO: U, V from vw, dw
return varname;
}
void CNRMIO::add_attributes_for_variable(const int& ncid, const int& varid, const std::string& varname)
{
if (varname == cf_latitude) {
if (varname == cf_time) { //HACK this should now be deprecated
ncpp::add_attribute(ncid, varid, "standard_name", CNRMIO::cf_time);
ncpp::add_attribute(ncid, varid, "long_name", CNRMIO::cf_time);
ncpp::add_attribute(ncid, varid, "units", "days since 1858-11-17 00:00:00");
} else if (varname == CNRMIO::cnrm_timestep) {
ncpp::add_attribute(ncid, varid, "long_name", "Forcing_Time_Step");
ncpp::add_attribute(ncid, varid, "units", "s");
} else if (varname == CNRMIO::cnrm_latitude) {
ncpp::add_attribute(ncid, varid, "standard_name", "latitude");
ncpp::add_attribute(ncid, varid, "long_name", "latitude");
ncpp::add_attribute(ncid, varid, "units", "degrees_north");
} else if (varname == cf_longitude) {
} else if (varname == CNRMIO::cnrm_longitude) {
ncpp::add_attribute(ncid, varid, "standard_name", "longitude");
ncpp::add_attribute(ncid, varid, "long_name", "longitude");
ncpp::add_attribute(ncid, varid, "units", "degrees_east");
} else if (varname == cf_altitude) {
} else if (varname == CNRMIO::cnrm_altitude) {
ncpp::add_attribute(ncid, varid, "standard_name", "altitude");
ncpp::add_attribute(ncid, varid, "long_name", "height above mean sea level");
ncpp::add_attribute(ncid, varid, "long_name", "altitude");
ncpp::add_attribute(ncid, varid, "units", "m");
ncpp::add_attribute(ncid, varid, "positive", "up");
ncpp::add_attribute(ncid, varid, "axis", "Z");
} else if (varname == cf_p) {
ncpp::add_attribute(ncid, varid, "standard_name", "air_pressure");
ncpp::add_attribute(ncid, varid, "long_name", "near surface air pressure");
ncpp::add_attribute(ncid, varid, "units", "Pa");
} else if (varname == cf_ta) {
ncpp::add_attribute(ncid, varid, "standard_name", "air_temperature");
ncpp::add_attribute