WSL/SLF GitLab Repository

Commit 3e406fd0 authored by Mathias Bavay's avatar Mathias Bavay
Browse files

Several bugs have been fixed in Coords: the setEPSG was not working, the...

Several bugs have been fixed in Coords: the setEPSG was not working, the altitude was not displayed with the << operator and the formatting of decimal_to_dms was broken in its handling of seconds...

The documentation for the resampling has been expanded and several other modules' documentation have been updated so their structure is built along the same lines. 

An error message in IOHandler (when failing to load a plugin) was not clear and has been changed.

Several bugs in SMEIO have been fixed (improper column names mapping, improper handling of easting/northing, unclear error message).

The example data set has been tagged as GMT+1 and the meteo_reading and 2D_interpolations modified accordingly (so the argument read on the command line is assumed to be in GMT+1 too)
parent 893b43d3
......@@ -20,7 +20,8 @@ int main(int /*argc*/, char** argv) {
DEMObject dem;
io.readDEM(dem);
//reading the meteorological data for the requested time step
//we assume that the time given on the command line is in TZ=+1
d1.setTimeZone(1.);
IOUtils::convertString(d1,argv[1]);
io.readMeteoData(d1, vecMeteo, vecStation);
......
......@@ -9,6 +9,7 @@ easting = 647900.000000
northing = 168780.000000
epsg = 21781
nodata = -999
tz = 1
fields = timestamp TA RH VW DW OSWR HS PSUM TSG TSS
[DATA]
2008-12-01T00:00 269.15 0.520 2.0 265 0 85.000 0.000 273.85 261.75
......
......@@ -9,6 +9,7 @@ easting = 791600.000000
northing = 180975.000000
epsg = 21781
nodata = -999
tz = 1
fields = timestamp TA RH VW DW OSWR HS TSG TSS
[DATA]
2008-12-01T00:00 265.35 1.000 6.0 134 0 36.000 273.15 265.75
......
......@@ -9,6 +9,7 @@ easting = 708900.000000
northing = 132850.000000
epsg = 21781
nodata = -999
tz = 1
fields = timestamp TA RH VW DW OSWR HS TSG TSS
[DATA]
2008-12-01T00:00 266.75 1.000 0.9 214 -1 158.000 273.45 266.85
......
......@@ -9,6 +9,7 @@ easting = 721610.000000
northing = 206300.000000
epsg = 21781
nodata = -999
tz = 1
fields = timestamp TA RH VW DW OSWR HS PSUM TSG TSS
[DATA]
2008-12-01T00:00 273.05 0.442 0.6 222 0 117.000 0.000 274.15 264.65
......
......@@ -9,6 +9,7 @@ easting = 552840.000000
northing = 115725.000000
epsg = 21781
nodata = -999
tz = 1
fields = timestamp TA RH VW DW OSWR HS PSUM TSG TSS
[DATA]
2008-12-01T00:00 269.45 0.369 1.2 171 0 31.000 0.000 273.55 263.35
......
......@@ -9,6 +9,7 @@ easting = 609450.000000
northing = 154250.000000
epsg = 21781
nodata = -999
tz = 1
fields = timestamp TA RH VW DW OSWR HS PSUM TSG TSS
[DATA]
2008-12-01T00:00 269.65 0.443 0.4 91 0 68.000 0.000 274.45 263.95
......
......@@ -9,6 +9,7 @@ easting = 699639.000000
northing = 167027.000000
epsg = 21781
nodata = -999
tz = 1
fields = timestamp TA RH VW DW OSWR HS PSUM TSG TSS
[DATA]
2008-12-01T00:00 265.95 0.988 0.0 335 0 142.000 0.000 273.55 264.45
......
......@@ -107,6 +107,16 @@ VW::arg1 = -0.5 80.
VW::filter2 = min
VW::arg2 = soft 0.
[Interpolations1D]
TA::resample = linear
RH::resample = linear
VW::resample = nearest_neighbour
VW::args = extrapolate
HNW::resample = linear
[Interpolations2D]
TA::algorithms = IDW_LAPSE CST_LAPSE
TA::cst_lapse = -0.008
......
......@@ -15,6 +15,8 @@ int main(int /*argc*/, char** argv) {
IOHandler raw_io(cfg);
BufferedIOHandler io(raw_io, cfg);
//we assume that the time given on the command line is in TZ=+1
d1.setTimeZone(1.);
IOUtils::convertString(d1,argv[1]);
io.readMeteoData(d1, vecMeteo, vecStation);
......
......@@ -138,6 +138,7 @@ Coords& Coords::operator=(const Coords& source) {
std::ostream& operator<<(std::ostream &os, const Coords& coord)
{
os << "<Coords>\n";
os << "Altitude\t" << coord.altitude << "\n";
os << "Lat/Long\t" << coord.printLatLon() << "\n";
os << "X/Y_coords\t" << "(" << coord.getEasting() << " , " << coord.getNorthing() << ")" << "\n";
os << "I/J_indices\t" << "(" << coord.getGridI() << " , " << coord.getGridJ() << ")" << "\n";
......@@ -566,37 +567,39 @@ short int Coords::getEPSG() const {
*/
void Coords::setEPSG(const short int epsg) {
bool found=false;
std::string coord_sys, coord_param;
if(!found && (epsg==21781)) {
coordsystem="CH1903";
coordparam="";
coord_sys="CH1903";
coord_param="";
found=true;
}
if(!found && (epsg>=32601) && (epsg<=32660)) {
//northern hemisphere
coordsystem="UTM";
coord_sys="UTM";
const short int zoneNumber = epsg-32600;
std::ostringstream osstream;
osstream << zoneNumber << "P";
coordparam=osstream.str();
coord_param=osstream.str();
found=true;
}
if(!found && (epsg>=32701) && (epsg<=32760)) {
//southern hemisphere
coordsystem="UTM";
coord_sys="UTM";
const short int zoneNumber = epsg-32700;
std::ostringstream osstream;
osstream << zoneNumber << "N";
coordparam=osstream.str();
coord_param=osstream.str();
found=true;
}
if(!found) {
//anything else has to be processed by proj4
coordsystem="PROJ4";
coord_sys="PROJ4";
std::ostringstream osstream;
osstream << epsg;
coordparam=osstream.str();
coord_param=osstream.str();
}
setProj(coord_sys, coord_param);
}
/////////////////////////////////////////////////////private methods
......@@ -720,9 +723,9 @@ std::string Coords::decimal_to_dms(const double& decimal) {
std::stringstream dms;
const int d = (int)floor(decimal);
const int m = (int)floor( (decimal - (double)d)*60. );
const int s = (int)floor( decimal - (double)d - (double)m/60. );
const double s = 3600.*(decimal - (double)d) - 60.*(double)m;
dms << d << "°" << m << "'" << s << "\"";
dms << d << "°" << m << "'" << fixed << setprecision(6) << s << "\"";
return dms.str();
}
......
......@@ -27,6 +27,26 @@ namespace mio {
* @section filters_modes Modes of operation
* It should be noted that filters often have two modes of operations: soft or hard. In soft mode, all value that is rejected is replaced by the filter parameter's value. This means that for a soft min filter set at 0.0, all values less than 0.0 will be replaced by 0.0. In hard mode, all rejected values are replaced by nodata.
*
* @section filtering_section Filtering section
* The filters are specified for each parameter in the [Filters] section. This section contains
* a list of the various meteo parameters with their associated choice of filtering algorithms and
* optional parameters.The filters are applied serialy, in the order they are given in. An example of such section is given below:
* @code
* [Filters]
* TA::filter1 = min_max
* TA::arg1 = 230 330
*
* RH::filter1 = min_max
* RH::arg1 = -0.2 1.2
* RH::filter2 = min_max
* RH::arg2 = soft 0.0 1.0
*
* HNW::filter1 = min
* HNW::arg1 = -0.1
* HNW::filter2 = min
* HNW::arg2 = soft 0.
* @endcode
*
* @section filters_available Available filters
* The filters that are currently available are the following:
* - rate: rate of change filter, see FilterAlgorithms::RateFilter
......
......@@ -183,7 +183,7 @@ IOInterface* IOHandler::getPlugin(const std::string& cfgkey, const std::string&
std::map<std::string, IOPlugin::IOPlugin>::iterator mapit = mapPlugins.find(op_src);
if (mapit == mapPlugins.end())
throw IOException(cfgkey + " does not seem to be valid plugin in file " + cfg.getSourceName(), AT);
throw IOException("Can not find plugin " + op_src + " as requested in file " + cfg.getSourceName() + ". Has its developer declared it in IOHandler::registerPlugins?", AT);
if ((mapit->second).io == NULL){
loadPlugin((mapit->second).libname, (mapit->second).classname, (mapit->second).dynLibrary, (mapit->second).io);
......
......@@ -65,6 +65,7 @@ namespace IOUtils {
//const double not_set = std::numeric_limits<double>::max()-2.;
const unsigned int unodata = (unsigned int)-1;
const int inodata = -999;
const short int snodata = -999;
const unsigned int npos = (unsigned int)-1; ///<npos is the out-of-range value
const double earth_radius = 6371e3; ///<Earth radius in meters
......
......@@ -42,28 +42,14 @@ class Meteo2DInterpolator; // forward declaration, cyclic header include
* memory of (eventual) previous time steps is kept. This means that all parameters and variables that are
* automatically calculated get recalculated anew for each time step.
*
* @section practical Practical use
* @section interpol2D_section Spatial interpolations section
* Practically, the user
* has to specify in his configuration file (typically io.ini), for each parameter to be interpolated, which
* spatial interpolations algorithms should be considered. This is provided as a space separated list of keywords
* spatial interpolations algorithms should be considered, in the [Interpolations2D] section. This is provided as a space separated list of keywords
* (one per interpolation algorithm). Please notice that some algorithms may require extra arguments.
* Then, each algorithm will be evaluated (through the use of its rating method) and receive a grade (that might
* depend on the number of available data, the quality of the data, etc). The algorithm that receives the higher
* score within the user list, will be used for interpolating the selected variable.
*
* @section keywords Available algorithms
* The keywords defining the algorithms are the following:
* - STD_PRESS: standard atmospheric pressure as a function of the elevation of each cell (see StandardPressureAlgorithm)
* - CST: constant value in each cell (see ConstAlgorithm)
* - CST_LAPSE: constant value reprojected to the elevation of the cell (see ConstLapseRateAlgorithm)
* - IDW: Inverse Distance Weighting averaging (see IDWAlgorithm)
* - IDW_LAPSE: Inverse Distance Weighting averaging with reprojection to the elevation of the cell (see IDWLapseAlgorithm)
* - RH: the dew point temperatures are interpolated using IDW_LAPSE, then reconverted locally to relative humidity (see RHAlgorithm)
* - WIND_CURV: the wind field (VW and DW) is interpolated using IDW_LAPSE and then altered depending on the local curvature and slope (taken from the DEM, see SimpleWindInterpolationAlgorithm)
* - USER: user provided grids to be read from disk (if available, see USERinterpolation) THIS IS NOT YET USABLE
*
* @section example Example of configuration file
* Here is an example of the interpolation section of an configuration file (io.ini):
* score within the user list, will be used for interpolating the selected variable. An example of such section is given below:
* @code
* [Interpolations2D]
* TA::algorithms = IDW_LAPSE CST_LAPSE
......@@ -78,6 +64,17 @@ class Meteo2DInterpolator; // forward declaration, cyclic header include
* P::algorithms = STD_PRESS
* @endcode
*
* @section keywords Available algorithms
* The keywords defining the algorithms are the following:
* - STD_PRESS: standard atmospheric pressure as a function of the elevation of each cell (see StandardPressureAlgorithm)
* - CST: constant value in each cell (see ConstAlgorithm)
* - CST_LAPSE: constant value reprojected to the elevation of the cell (see ConstLapseRateAlgorithm)
* - IDW: Inverse Distance Weighting averaging (see IDWAlgorithm)
* - IDW_LAPSE: Inverse Distance Weighting averaging with reprojection to the elevation of the cell (see IDWLapseAlgorithm)
* - RH: the dew point temperatures are interpolated using IDW_LAPSE, then reconverted locally to relative humidity (see RHAlgorithm)
* - WIND_CURV: the wind field (VW and DW) is interpolated using IDW_LAPSE and then altered depending on the local curvature and slope (taken from the DEM, see SimpleWindInterpolationAlgorithm)
* - USER: user provided grids to be read from disk (if available, see USERinterpolation) THIS IS NOT YET USABLE
*
* @section lapse Lapse rates
* Several algorithms use elevation trends, currently modelled as a linear relation. The slope of this linear relation can
* sometimes be provided by the end user (through his io.ini configuration file), otherwise it is computed from the data.
......@@ -107,6 +104,7 @@ class Meteo2DInterpolator; // forward declaration, cyclic header include
* The interpolation algorithms have been inspired by the following papers:
* - "A Meteorological Distribution System for High-Resolution Terrestrial Modeling (MicroMet)", Liston and Elder, Journal of Hydrometeorology 7 (2006), 217-234.
* - "Simulating wind fields and snow redistribution using terrain-based parameters to model snow accumulation and melt over a semi-arid mountain catchment", Adam Winstral and Danny Marks, Hydrol. Process. 16 (2002), 3585– 3603. DOI: 10.1002/hyp.1238 [NOT YET IMPLEMENTED]
* - "Geostatistics for Natural Resources Evaluation", Pierre Goovaerts, Oxford University Press, Applied Geostatistics Series, 1997, 483 p., ISBN 0-19-511538-4
*
* @author Mathias Bavay
* @date 2010-04-12
......
......@@ -25,6 +25,23 @@ namespace mio {
* The resampling infrastructure is described in ResamplingAlgorithms (for its API).
* The goal of this page is to give an overview of the available resampling algorithms and their usage.
*
* @section resampling_section Resampling section
* The resampling is specified for each parameter in the [Interpol1D] section. This section contains
* a list of the various meteo parameters with their associated choice of resampling algorithm and
* optional parameters. If a meteo parameter is not listed in this section, a linear resampling would be
* assumed. An example of such section is given below:
* @code
* [Interpolations1D]
* TA::resample = linear
*
* RH::resample = linear
*
* VW::resample = nearest_neighbour
* VW::args = extrapolate
*
* HNW::resample = linear
* @endcode
*
* @section algorithms_available Available Resampling Algorithms
* Two algorithms for the resampling are implemented:
* - linear: linear data resampling, see ResamplingAlgorithms::LinearResampling
......
......@@ -31,7 +31,7 @@ PROJECT_NAME = MeteoIODoc
# This could be handy for archiving the generated documentation or
# if some version control system is used.
PROJECT_NUMBER = MeteoIODoc-1.0
PROJECT_NUMBER = MeteoIODoc-1.1
# The OUTPUT_DIRECTORY tag is used to specify the (relative or absolute)
# base path where the generated documentation will be put.
......
......@@ -60,7 +60,16 @@ bool SMETIO::initStaticData()
}
double& SMETIO::getParameter(const std::string& columnName, MeteoData& md)
{
{//HACK: the whole name mapping is a big hack. Replace with proper mapping!!!
if(columnName=="OSWR") {
MeteoData::Parameters paramindex = mapParameterByName["RSWR"];
return md.param(paramindex);
}
if(columnName=="PSUM") {
MeteoData::Parameters paramindex = mapParameterByName["RSWR"];
return md.param(paramindex);
}
MeteoData::Parameters paramindex = mapParameterByName[columnName];
return md.param(paramindex);
}
......@@ -469,9 +478,10 @@ void SMETIO::readHeader(const char& eoln, const std::string& filename, bool& loc
IOUtils::stripComments(line);
IOUtils::trim(line);
if (line != "")
if (line != "") {
if (!IOUtils::readKeyValuePair(line, "=", mapHeader))
throw InvalidFormatException("Invalid key value pair in section [Header]", AT);
}
}
//Now extract info from mapHeader
......@@ -479,15 +489,29 @@ void SMETIO::readHeader(const char& eoln, const std::string& filename, bool& loc
IOUtils::getValueForKey(mapHeader, "station_name", sd.stationName, IOUtils::nothrow);
IOUtils::getValueForKey(mapHeader, "tz", timezone, IOUtils::nothrow);
double lat=IOUtils::nodata, lon=IOUtils::nodata, alt=IOUtils::nodata;
//trying to read easting/northing
double easting=IOUtils::nodata, northing=IOUtils::nodata, alt=IOUtils::nodata;
short int epsg=IOUtils::snodata;
IOUtils::getValueForKey(mapHeader, "easting", easting, IOUtils::nothrow);
if (easting != IOUtils::nodata){ //HACK
IOUtils::getValueForKey(mapHeader, "northing", northing);
IOUtils::getValueForKey(mapHeader, "altitude", alt);
IOUtils::getValueForKey(mapHeader, "epsg", epsg);
sd.position.setEPSG(epsg);
sd.position.setXY(easting, northing, alt);
locationInHeader = true;
} else {
locationInHeader = false;
}
//now trying to read lat/long (second, so that it has precedence over east/north coordinates)
double lat=IOUtils::nodata, lon=IOUtils::nodata;
IOUtils::getValueForKey(mapHeader, "latitude", lat, IOUtils::nothrow);
if (lat != IOUtils::nodata){ //HACK
IOUtils::getValueForKey(mapHeader, "longitude", lon);
IOUtils::getValueForKey(mapHeader, "altitude", alt);
sd.position.setLatLon(lat, lon, alt);
locationInHeader = true;
} else {
locationInHeader = false;
}
IOUtils::getValueForKey(mapHeader, "fields", vecDataSequence);
......@@ -519,7 +543,7 @@ void SMETIO::checkSignature(const std::vector<std::string>& vecSignature, const
else if (type == "BINARY")
isAscii = false;
else
throw InvalidFormatException("The 3rd column in the file " + filename + " must be either ASCII or BINARY", AT);
throw InvalidFormatException("The 3rd column in the signature of file " + filename + " must be either ASCII or BINARY", AT);
}
void SMETIO::writeMeteoData(const std::vector< std::vector<MeteoData> >& vecMeteo,
......
Supports Markdown
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