WSL/SLF GitLab Repository

ProcessingBlock.cc 9.11 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
/***********************************************************************************/
/*  Copyright 2009 WSL Institute for Snow and Avalanche Research    SLF-DAVOS      */
/***********************************************************************************/
/* This file is part of MeteoIO.
    MeteoIO is free software: you can redistribute it and/or modify
    it under the terms of the GNU Lesser General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    MeteoIO is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU Lesser General Public License for more details.

    You should have received a copy of the GNU Lesser General Public License
    along with MeteoIO.  If not, see <http://www.gnu.org/licenses/>.
*/
18
#include <meteoio/meteofilters/ProcessingBlock.h>
19
20
#include <meteoio/meteofilters/FilterMin.h>
#include <meteoio/meteofilters/FilterMax.h>
21
22
23
#include <meteoio/meteofilters/FilterMinMax.h>
#include <meteoio/meteofilters/FilterMeanAvg.h>
#include <meteoio/meteofilters/FilterMedianAvg.h>
24
#include <meteoio/meteofilters/FilterWindAvg.h>
25
#include <meteoio/meteofilters/FilterStdDev.h>
26
#include <meteoio/meteofilters/FilterRate.h>
27
#include <meteoio/meteofilters/FilterUnheatedHNW.h>
28
#include <meteoio/meteofilters/FilterTukey.h>
29
#include <meteoio/meteofilters/FilterMAD.h>
30
#include <meteoio/meteofilters/ProcUndercatch_WMO.h>
31
#include <meteoio/meteofilters/ProcUnventilatedT.h>
32
33
#include <meteoio/meteofilters/ProcAdd.h>
#include <meteoio/meteofilters/ProcMult.h>
34
35
#include <meteoio/meteofilters/ProcExpSmoothing.h>
#include <meteoio/meteofilters/ProcWMASmoothing.h>
36
37

namespace mio {
38
39
40
41
42
43
44
/**
 * @page processing Processing overview
 * The pre-processing infrastructure is described in ProcessingBlock (for its API). The goal of this page is to give an overview of the available filters and processing elements and their usage.
 *
 * @section processing_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.
 *
45
 * Several filter take arguments describing a processing window (for example, FilterStdDev). In such a case, two values are given:
46
 * - the minimum time span of the window
47
48
 * - the minimum number of points in the window
 *
49
50
 * The ProcessingBlock will walk through the data, starting at the current point and adding points to the processing window.  When both of these criterias are met,
 * the window is accepted. This means that a window defined as "6 21600" is a window that contains 6 points minimum AND spans at least 21600 seconds.
51
 *
52
53
54
55
56
57
58
59
 * @section processing_section Filtering section
 * The filters are specified for each parameter in the [Filters] section. This section contains
 * a list of the various meteo parameters (see MeteoData) 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
Mathias Bavay's avatar
Mathias Bavay committed
60
 *
61
62
63
64
 * RH::filter1	= min_max
 * RH::arg1	= -0.2 1.2
 * RH::filter2	= min_max
 * RH::arg2	= soft 0.0 1.0
Mathias Bavay's avatar
Mathias Bavay committed
65
 *
66
67
68
69
70
71
72
 * HNW::filter1	= min
 * HNW::arg1	= -0.1
 * HNW::filter2	= min
 * HNW::arg2	= soft 0.
 * @endcode
 *
 * @section processing_available Available processing elements
73
 * New filters can easily be developed. The filters that are currently available are the following:
74
 * - RATE: rate of change filter, see FilterRate
75
76
77
 * - MIN_MAX: range check filter, see FilterMinMax
 * - MIN: minimum check filter, see FilterMin
 * - MAX: maximum check filter, see FilterMax
Mathias Bavay's avatar
Mathias Bavay committed
78
 * - STD_DEV: reject data outside mean +/- k*stddev, see FilterStdDev
79
 * - MAD: median absolute deviation, see FilterMAD
80
 * - TUKEY: Tukey53H spike detection, based on median, see FilterTukey
81
 * - UNHEATED_RAINGAUGE: detection of snow melting in a rain gauge, see FilterUnheatedHNW
82
83
 *
 * A few data transformations are also supported besides filtering:
84
 * - EXP_SMOOTHING: exponential smoothing of data, see ProcExpSmoothing
Mathias Bavay's avatar
Mathias Bavay committed
85
 * - WMA_SMOOTHING weighted moving average smoothing of data, see ProcWMASmoothing
86
87
 * - MEDIAN_AVG: running median average over a given window, see FilterMedianAvg
 * - MEAN_AVG: running mean average over a given window, see FilterMeanAvg
88
 * - WIND_AVG: vector average over a given window, see FilterWindAvg (currently, getting both vw AND dw is broken)
89
90
 * - ADD: adds a given offset to the data, see ProcAdd
 * - MULT: multiply the data by a given factor, see ProcMult
91
 * - UNDERCATCH_WMO: WMO rain gauge correction for undercatch, using various correction models, see ProcUndercatch_WMO
92
 * - UNVENTILATED_T: unventilated temperature sensor correction, see ProcUnventilatedT
93
 *
94
 */
95
96
97
98
99
100

std::set<std::string> BlockFactory::availableBlocks;
const bool BlockFactory::__init = BlockFactory::initStaticData();

bool BlockFactory::initStaticData()
{
101
102
	availableBlocks.insert("MIN");
	availableBlocks.insert("MAX");
103
104
	availableBlocks.insert("MIN_MAX");
	availableBlocks.insert("MEAN_AVG");
105
	availableBlocks.insert("MEDIAN_AVG");
106
	availableBlocks.insert("WIND_AVG");
107
	availableBlocks.insert("STD_DEV");
108
	availableBlocks.insert("RATE");
109
	availableBlocks.insert("TUKEY");
110
	availableBlocks.insert("MAD");
111
	availableBlocks.insert("UNHEATED_RAINGAUGE");
112
	availableBlocks.insert("UNDERCATCH_WMO");
113
	availableBlocks.insert("UNVENTILATED_T");
114
115
	availableBlocks.insert("ADD");
	availableBlocks.insert("MULT");
116
117
	availableBlocks.insert("EXP_SMOOTHING");
	availableBlocks.insert("WMA_SMOOTHING");
118
119
120
121
122
123
124
125
126
	return true;
}

ProcessingBlock* BlockFactory::getBlock(const std::string& blockname, const std::vector<std::string>& vec_args)
{
	//Check whether algorithm theoretically exists
	if (availableBlocks.find(blockname) == availableBlocks.end())
		throw UnknownValueException("The processing block '"+blockname+"' does not exist" , AT);

127
128
129
130
131
	if (blockname == "MIN"){
		return new FilterMin(vec_args);
	} else if (blockname == "MAX"){
		return new FilterMax(vec_args);
	} else if (blockname == "MIN_MAX"){
132
133
134
		return new FilterMinMax(vec_args);
	} else if (blockname == "MEAN_AVG"){
		return new FilterMeanAvg(vec_args);
135
136
	} else if (blockname == "MEDIAN_AVG"){
		return new FilterMedianAvg(vec_args);
137
138
	} else if (blockname == "WIND_AVG"){
		return new FilterWindAvg(vec_args);
139
140
	} else if (blockname == "STD_DEV"){
		return new FilterStdDev(vec_args);
141
	} else if (blockname == "RATE"){
142
		return new FilterRate(vec_args);
143
144
	} else if (blockname == "TUKEY"){
		return new FilterTukey(vec_args);
145
146
	} else if (blockname == "MAD"){
		return new FilterMAD(vec_args);
147
148
	} else if (blockname == "UNHEATED_RAINGAUGE"){
		return new FilterUnheatedHNW(vec_args);
149
150
	} else if (blockname == "UNDERCATCH_WMO"){
		return new ProcUndercatch_WMO(vec_args);
151
152
	} else if (blockname == "UNVENTILATED_T"){
		return new ProcUnventilatedT(vec_args);
153
154
155
156
	} else if (blockname == "MULT"){
		return new ProcMult(vec_args);
	} else if (blockname == "ADD"){
		return new ProcAdd(vec_args);
157
158
159
160
	} else if (blockname == "EXP_SMOOTHING"){
		return new ProcExpSmoothing(vec_args);
	} else if (blockname == "WMA_SMOOTHING"){
		return new ProcWMASmoothing(vec_args);
161
	} else {
Mathias Bavay's avatar
Mathias Bavay committed
162
		throw IOException("The processing block '"+blockname+"' has not been declared! " , AT);
163
164
165
166
167
	}

}

ProcessingBlock::ProcessingBlock(const std::string& name) : block_name(name)
168
{}
169

170
171
172
173
174
175
void ProcessingBlock::convert_args(const unsigned int& min_nargs, const unsigned int& max_nargs,
                               const std::vector<std::string>& vec_args, std::vector<double>& dbl_args)
{
	if ((vec_args.size() < min_nargs) || (vec_args.size() > max_nargs))
		throw InvalidArgumentException("Wrong number of arguments for filter/processing element \"" + getName() + "\"", AT);

176
	for (size_t ii=0; ii<vec_args.size(); ii++){
177
178
179
180
181
182
		double tmp = IOUtils::nodata;
		IOUtils::convertString(tmp, vec_args[ii]);
		dbl_args.push_back(tmp);
	}
}

183
std::string ProcessingBlock::getName() const {
184
185
186
187
188
	return block_name;
}

ProcessingBlock::~ProcessingBlock() {}

189
std::ostream& operator<<(std::ostream& os, const ProcessingBlock& data) {
190
	os << "[" << data.block_name << " ";
191
	os << data.properties;
192
193
194
195
	os << "]";
	return os;
}

196
const ProcessingProperties& ProcessingBlock::getProperties() const {
197
198
199
	return properties;
}

200
std::ostream& operator<<(std::ostream& os, const ProcessingProperties& data) {
201
202
203
204
205
206
207
208
	const double h_before = data.time_before.getJulianDate()*24.;
	const double h_after = data.time_after.getJulianDate()*24.;
	const unsigned int p_before = data.points_before;
	const unsigned int p_after = data.points_after;

	os << "{";
	if(h_before>0. || h_after>0.) os << "-" << h_before << " +" << h_after << " h; ";
	if(p_before>0 || p_after>0) os << "-" << p_before << " +" << p_after << " pts; ";
209
210
211
	if(data.stage==ProcessingProperties::first)
		os << "p¹";
	if(data.stage==ProcessingProperties::second)
212
		os << "p²";
213
214
	if(data.stage==ProcessingProperties::both)
		os << "p½";
215
216
217
218
	os << "}";
	return os;
}

219
} //end namespace