WSL/SLF GitLab Repository

Commit 10e1caff authored by Mathias Bavay's avatar Mathias Bavay
Browse files

Now multiple stations are properly supported by the SHADE filter (when...

Now multiple stations are properly supported by the SHADE filter (when computing the mask from the DEM)
parent 328a1236
......@@ -41,7 +41,7 @@ const double ProcShade::snow_albedo = .85; //snow
const double ProcShade::snow_thresh = .1; //if snow height greater than this threshold -> snow albedo
ProcShade::ProcShade(const std::vector<std::string>& vec_args, const std::string& name, const Config &i_cfg)
: ProcessingBlock(name), cfg(i_cfg), dem(), Suns(), mask()
: ProcessingBlock(name), cfg(i_cfg), dem(), Suns(), masks()
{
parse_args(vec_args);
properties.stage = ProcessingProperties::first; //for the rest: default values
......@@ -53,9 +53,10 @@ void ProcShade::process(const unsigned int& param, const std::vector<MeteoData>&
ovec = ivec;
if (ovec.size()==0) return;
//check if the station already has an associated SunObject
const string stationHash = ovec[0].meta.getHash();
map< string, SunObject >::iterator it = Suns.find( stationHash );
//check if the station already has an associated SunObject
std::map< std::string, SunObject >::iterator it = Suns.find( stationHash );
if (it==Suns.end()) {
const Coords position( ovec[0].meta.position );
const SunObject tmp(position.getLat(), position.getLon(), position.getAltitude());
......@@ -63,12 +64,20 @@ void ProcShade::process(const unsigned int& param, const std::vector<MeteoData>&
it = Suns.find( stationHash );
}
//compute the mask from DEM if a DEM is available
if (mask.empty()) {
Coords position( ovec[0].meta.position );
if (!dem.gridify(position))
throw NoDataException("In filter '"+block_name+"', station '"+ovec[0].meta.stationID+"' is not included in the DEM", AT);
dem.getHorizon(position, 10., mask);
//check if the station already has an associated mask, first as wildcard then by station hash
std::map< std::string , std::vector< std::pair<double,double> > >::iterator mask = masks.find( "*" );
if (mask==masks.end()) {
//now look for our specific station hash
mask = masks.find( stationHash );
if (mask==masks.end()) {
Coords position( ovec[0].meta.position );
if (!dem.gridify(position))
throw NoDataException("In filter '"+block_name+"', station '"+ovec[0].meta.stationID+"' is not included in the DEM", AT);
std::vector< std::pair<double,double> > tmp_mask;
dem.getHorizon(position, 10., tmp_mask);
masks[ stationHash ] = tmp_mask;
mask = masks.find( stationHash);
}
}
for (size_t ii=0; ii<ovec.size(); ii++){
......@@ -79,7 +88,7 @@ void ProcShade::process(const unsigned int& param, const std::vector<MeteoData>&
double sun_azi, sun_elev;
it->second.position.getHorizontalCoordinates(sun_azi, sun_elev);
const double mask_elev = getMaskElevation(sun_azi);
const double mask_elev = getMaskElevation(mask->second, sun_azi);
if (mask_elev>0 && mask_elev>sun_elev) { //in the shade
const double TA=ovec[ii](MeteoData::TA), RH=ovec[ii](MeteoData::RH), HS=ovec[ii](MeteoData::HS), RSWR=ovec[ii](MeteoData::RSWR);
double ISWR=ovec[ii](MeteoData::ISWR);
......@@ -111,7 +120,7 @@ void ProcShade::process(const unsigned int& param, const std::vector<MeteoData>&
}
//linear interpolation between the available points
double ProcShade::getMaskElevation(const double& azimuth) const
double ProcShade::getMaskElevation(const std::vector< std::pair<double,double> > &mask, const double& azimuth) const
{
const std::vector< std::pair<double, double> >::const_iterator next = std::upper_bound(mask.begin(), mask.end(), make_pair(azimuth, 0.), sort_pred()); //first element that is > azimuth
......@@ -202,7 +211,9 @@ void ProcShade::parse_args(const std::vector<std::string>& vec_args)
const std::string prefix = ( IOUtils::isAbsolutePath(in_filename) )? "" : root_path+"/";
const std::string path = IOUtils::getPath(prefix+in_filename, true); //clean & resolve path
const std::string filename = path + "/" + IOUtils::getFilename(in_filename);
std::vector< std::pair<double,double> > mask;
readMask(getName(), filename, mask);
masks["*"] = mask; //this mask is valid for ALL stations
} else
throw InvalidArgumentException("Wrong number of arguments for filter " + getName(), AT);
......
......@@ -30,11 +30,13 @@ namespace mio {
/**
* @class ProcShade
* @ingroup processing
* @author Mathias Bavay
* @date 2016-01-22
* @brief Apply a shading mask to the Incoming or Reflected Short Wave Radiation
* A shading mask (that will be linearly interpolated between the provided points) must be provided in a separate file,
* simply containing the horizon elevation (in deg.) as a function of azimuth (in deg.):
* A shading mask that is either computed from the DEM or read from a separate file will be applied to the radiation
* and combined with the radiation splitting model in order to properly compute the shading effects on the measurement point.
* This mask will be linearly interpolated between the provided points in order to be applied to the true sun position.
*
* When providing thge shading mask in a separate file, the same mask will be applied to all stations. It simply need to
* contain the horizon elevation (in deg.) as a function of azimuth (in deg.):
* @code
* 0 5
* 15 25
......@@ -50,8 +52,8 @@ namespace mio {
* @endcode
*
* If no arguments are provided, then it will compute the mask from the Digital Elevation Model. In such as case,
* a DEM must be declared in the [Input] section and must contain the stations of interest. Please make sure that the extend of the
* DEM is appropriate to correctly compute the shading effects!
* a DEM must be declared in the [Input] section and must contain the stations of interest as a mask will be computed for
* each station. Please make sure that the extend of the DEM is appropriate to correctly compute the shading effects!
*/
class ProcShade : public ProcessingBlock {
......@@ -65,12 +67,12 @@ class ProcShade : public ProcessingBlock {
static void readMask(const std::string& filter, const std::string& filename, std::vector< std::pair<double,double> > &o_mask);
void parse_args(const std::vector<std::string>& vec_args);
double getMaskElevation(const double& azimuth) const;
double getMaskElevation(const std::vector< std::pair<double,double> > &mask, const double& azimuth) const;
const Config &cfg;
DEMObject dem;
std::map<std::string, SunObject> Suns;
std::vector< std::pair<double,double> > mask;
std::map< std::string , std::vector< std::pair<double,double> > > masks;
static const double soil_albedo, snow_albedo, snow_thresh; ///< parametrize the albedo from HS
};
......
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