WindowedFilter.cc

/***********************************************************************************/
/*  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/>.
*/
#include <meteoio/meteofilters/WindowedFilter.h>

using namespace std;

namespace mio {

WindowedFilter::WindowedFilter(const std::string& name)
	: FilterBlock(name), min_time_span(0.0, 0.), centering(WindowedFilter::center),
	  last_start(0), last_end(0), min_data_points(1), vec_window(), is_soft(false)
{}

unsigned int WindowedFilter::get_centering(std::vector<std::string>& vec_args)
{
	if (!vec_args.empty()){
		if (vec_args.front() == "left"){
			vec_args.erase(vec_args.begin());
			return WindowedFilter::left;
		} else if (vec_args.front() == "right"){
			vec_args.erase(vec_args.begin());
			return WindowedFilter::right;
		} else if (vec_args.front() == "center"){
			vec_args.erase(vec_args.begin());
			return WindowedFilter::center;
		}
	}

	return WindowedFilter::center; //the default
}

/**
 * @brief A function that cuts out the desired window for the 'index' element within
 *        ivec, the window elements are stored into vec_window
 *        Calls to this function have to start with index 0, then 1, 2, 3, ...
 *        vec_window is not allowed to be changed between two calls
 * @param index The index of the element in ivec that requires a window
 * @param ivec The original sequence of data points
 */
const std::vector<const MeteoData*>& WindowedFilter::get_window(const size_t& index,
                                                                const std::vector<MeteoData>& ivec)
{
	if(index==0) {
		vec_window.clear();
		vec_window.reserve(min_data_points*2); //to have enough margin for the time criteria
	}

	size_t start, end;
	if(!get_window_specs(index, ivec, start, end)) {
		vec_window.clear();
		vec_window.reserve(min_data_points*2); //to have enough margin for the time criteria
		vec_window.push_back( &ivec[index] );
	}

	if(index==0) {
		for(size_t ii=start; ii<=end; ii++) vec_window.push_back( &ivec[ii] );
	} else {
		if(last_start<start) {
			vec_window.erase( vec_window.begin(),  vec_window.begin()+(start-last_start));
		}
		if(last_start>start) {
			const std::vector<const MeteoData*>::iterator it = vec_window.begin();
			for(size_t ii=start; ii<=last_start; ii++) vec_window.insert(it, &ivec[ii]);
		}
		if(last_end<end) {
			for(size_t ii=last_end+1; ii<=end; ii++) vec_window.push_back( &ivec[ii] );
		}
		if(last_end>end) {
			vec_window.erase( vec_window.end()-(last_end-end), vec_window.end() );
		}
	}

	//save window metadata
	last_start = start;
	last_end = end;
	return vec_window;
}

/**
 * @brief A function that computes the start and end for a window for the 'index' element from ivec
 * @param index The index of the element in ivec that requires a window
 * @param ivec The original sequence of data points
 * @param start the start index of the window
 * @param end the end index of the window
 * @return true if success, false if a window could not be computed
 */
bool WindowedFilter::get_window_specs(const size_t& index, const std::vector<MeteoData>& ivec, size_t& start, size_t& end) const
{	/*
	The principle is too compute the first index that matches the minimum number of points criteria,
	the the one that matches the minimum time window,
	then combine them (with the equivalent of OR: we take the MIN index).
	Afterward, we compute the last index [...] for number of points
	and the last index [...] for the time window
	and combine them (with the equivalent of OR: we take the MIN index).
	(or vice versa for right centering)
	*/
	const Date date = ivec[index].date;
	start = end = index; //for proper initial value, so we can bail out without worries

	if(centering == WindowedFilter::left) {
		//get start of window
		size_t start_elements = min_data_points - 1; //start elements criteria
		if(start_elements>index) {
			if(!is_soft) return false;
			start_elements = index; //as many as possible
		}
		const Date start_date = date - min_time_span;
		size_t start_time_idx = IOUtils::seek(start_date, ivec, false); //start time criteria
		if(start_time_idx!=IOUtils::npos) start_time_idx = (start_time_idx>0)?start_time_idx-1:IOUtils::npos;
		if(start_time_idx==IOUtils::npos) {
			if(!is_soft) return false;
			start_time_idx=0; //first possible element
		}
		const size_t elements_left = MAX(index - start_time_idx, start_elements);
		start = index - elements_left;

		//get end of window
		if(!is_soft) return true; //with end=index
		size_t end_elements = (min_data_points>(elements_left+1))?min_data_points - (elements_left + 1):0;
		const Date end_date = ivec[start].date+min_time_span;
		size_t end_time_idx = (end_date>date)?IOUtils::seek(end_date, ivec, false):index; //end time criteria
		if(end_time_idx==IOUtils::npos) {
			if(!is_soft) return false;
			end_time_idx=ivec.size()-1; //last possible element
		}
		const size_t elements_right = MAX(end_time_idx - index, end_elements);
		end = index + elements_right;
	}

	if(centering == WindowedFilter::right) {
		//get end of window
		size_t end_elements = min_data_points - 1; //end elements criteria
		if(end_elements>(ivec.size()-1-index)) {
			if(!is_soft) return false;
			end_elements = (ivec.size()-1-index); //as many as possible
		}
		const Date end_date = date+min_time_span;
		size_t end_time_idx = IOUtils::seek(end_date, ivec, false); //end time criteria
		if(end_time_idx==IOUtils::npos) {
			if(!is_soft) return false;
			end_time_idx=ivec.size()-1; //last possible element
		}
		const size_t elements_right = MAX(end_time_idx - index, end_elements);
		end = index + elements_right;

		//get start of window
		if(!is_soft) return true; //with start=index
		size_t start_elements = (min_data_points>(elements_right+1))?min_data_points - (elements_right + 1):0;
		const Date start_date = ivec[end].date-min_time_span;
		size_t start_time_idx = (start_date<date)?IOUtils::seek(start_date, ivec, false):index; //start time criteria
		if(start_time_idx!=IOUtils::npos && start_time_idx!=index) start_time_idx = (start_time_idx>0)?start_time_idx-1:IOUtils::npos;
		if(start_time_idx==IOUtils::npos) {
			if(!is_soft) return false;
			end_time_idx=0; //first possible element
		}
		const size_t elements_left = MAX(index - start_time_idx, start_elements);
		start = index - elements_left;
	}

	if(centering == WindowedFilter::center) {
		//get start of ideal window
		size_t start_elements = min_data_points/2; //start elements criteria
		if(start_elements>index) {
			if(!is_soft) return false;
			start_elements = index; //as many as possible
		}
		const Date start_date = date - min_time_span/2;
		size_t start_time_idx = IOUtils::seek(start_date, ivec, false); //start time criteria
		if(start_time_idx!=IOUtils::npos) start_time_idx = (start_time_idx>0)?start_time_idx-1:IOUtils::npos;
		if(start_time_idx==IOUtils::npos) {
			if(!is_soft) return false;
			start_time_idx=0; //first possible element
		}
		const size_t elements_left = MAX(index - start_time_idx, start_elements);
		start = index - elements_left;

		//get end of ideal window
		size_t end_elements = min_data_points/2; //end elements criteria
		if(end_elements>(ivec.size()-1-index)) {
			if(!is_soft) return false;
			end_elements = (ivec.size()-1-index); //as many as possible
		}
		const Date end_date = date+min_time_span/2;
		size_t end_time_idx = IOUtils::seek(end_date, ivec, false); //end time criteria
		if(end_time_idx==IOUtils::npos) {
			if(!is_soft) return false;
			end_time_idx=ivec.size()-1; //last possible element
		}
		const size_t elements_right = MAX(end_time_idx - index, end_elements);
		end = index + elements_right;

		//now, check (and modify) if the window could not be centered
		if(elements_left==elements_right) return true;
		if(!is_soft) return false;
		if(elements_left<elements_right) { //we hit the left border
			//get again the end of window
			size_t end_elems = (min_data_points>(elements_left+1))?min_data_points - (elements_left + 1):0;
			const Date end_dt = ivec[start].date+min_time_span;
			size_t end_tm_idx = (end_dt>date)?IOUtils::seek(end_dt, ivec, false):index; //end time criteria
			if(end_tm_idx==IOUtils::npos) {
				end_tm_idx=ivec.size()-1; //last possible element
			}
			const size_t elems_right = MAX(end_tm_idx - index, end_elems);
			end = index + elems_right;
		} else { //we hit the right border
			//get again the start of window
			size_t start_elems = (min_data_points>(elements_right+1))?min_data_points - (elements_right + 1):0;
			const Date start_dt = ivec[end].date-min_time_span;
			size_t start_tm_idx = (start_dt<date)?IOUtils::seek(start_dt, ivec, false):index; //start time criteria
			if(start_tm_idx!=IOUtils::npos && start_tm_idx!=index) start_tm_idx = (start_tm_idx>0)?start_tm_idx-1:IOUtils::npos;
			if(start_tm_idx==IOUtils::npos) {
				end_time_idx=0; //first possible element
			}
			const size_t elems_left = MAX(index - start_tm_idx, start_elems);
			start = index - elems_left;
		}
	}

	return true;
}

} //namespace