pktools 2.6.7
Processing Kernel for geospatial data
Public Types | Public Member Functions | Static Public Member Functions | List of all members
statfactory::StatFactory Class Reference

Public Types

enum  INTERPOLATION_TYPE {
  undefined =0 , linear =1 , polynomial =2 , cspline =3 ,
  cspline_periodic =4 , akima =5 , akima_periodic =6
}
 
enum  DISTRIBUTION_TYPE { uniform =1 , gaussian =2 }
 

Public Member Functions

INTERPOLATION_TYPE getInterpolationType (const std::string interpolationType)
 
DISTRIBUTION_TYPE getDistributionType (const std::string distributionType)
 
void getNodataValues (std::vector< double > &nodatav) const
 
bool isNoData (double value) const
 
unsigned int pushNodDataValue (double noDataValue)
 
unsigned int setNoDataValues (std::vector< double > vnodata)
 
double getRandomValue (const gsl_rng *r, const std::string type, double a=0, double b=1) const
 
template<class T >
mymin (const typename std::vector< T > &v) const
 
template<class T >
mymax (const typename std::vector< T > &v) const
 
template<class T >
mymin (const typename std::vector< T > &v, T minConstraint) const
 
template<class T >
mymax (const typename std::vector< T > &v, T maxConstraint) const
 
template<class T >
absmin (const typename std::vector< T > &v) const
 
template<class T >
absmax (const typename std::vector< T > &v) const
 
template<class T >
std::vector< T >::const_iterator mymin (const typename std::vector< T > &v, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end) const
 
template<class T >
std::vector< T >::iterator mymin (const typename std::vector< T > &v, typename std::vector< T >::iterator begin, typename std::vector< T >::iterator end) const
 
template<class T >
std::vector< T >::const_iterator mymin (const typename std::vector< T > &v, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end, T minConstraint) const
 
template<class T >
std::vector< T >::iterator mymin (const typename std::vector< T > &v, typename std::vector< T >::iterator begin, typename std::vector< T >::iterator end, T minConstraint) const
 
template<class T >
std::vector< T >::const_iterator mymax (const std::vector< T > &v, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end) const
 
template<class T >
std::vector< T >::iterator mymax (const std::vector< T > &v, typename std::vector< T >::iterator begin, typename std::vector< T >::iterator end) const
 
template<class T >
std::vector< T >::const_iterator mymax (const std::vector< T > &v, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end, T maxConstraint) const
 
template<class T >
std::vector< T >::iterator mymax (const std::vector< T > &v, typename std::vector< T >::iterator begin, typename std::vector< T >::iterator end, T maxConstraint) const
 
template<class T >
std::vector< T >::const_iterator absmin (const std::vector< T > &v, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end) const
 
template<class T >
std::vector< T >::const_iterator absmax (const std::vector< T > &v, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end) const
 
template<class T >
void minmax (const std::vector< T > &v, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end, T &theMin, T &theMax) const
 
template<class T >
sum (const std::vector< T > &v) const
 
template<class T >
double mean (const std::vector< T > &v) const
 
template<class T >
void eraseNoData (std::vector< T > &v) const
 
template<class T >
unsigned int nvalid (const std::vector< T > &v) const
 
template<class T >
median (const std::vector< T > &v) const
 
template<class T >
double var (const std::vector< T > &v) const
 
template<class T >
double moment (const std::vector< T > &v, int n) const
 
template<class T >
double cmoment (const std::vector< T > &v, int n) const
 
template<class T >
double skewness (const std::vector< T > &v) const
 
template<class T >
double kurtosis (const std::vector< T > &v) const
 
template<class T >
void meanVar (const std::vector< T > &v, double &m1, double &v1) const
 
template<class T1 , class T2 >
void scale2byte (const std::vector< T1 > &input, std::vector< T2 > &output, unsigned char lbound=0, unsigned char ubound=255) const
 
template<class T >
void distribution (const std::vector< T > &input, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end, std::vector< double > &output, int nbin, T &minimum, T &maximum, double sigma=0, const std::string &filename="") const
 
template<class T >
void distribution (const std::vector< T > &input, std::vector< double > &output, int nbin, double sigma=0, const std::string &filename="") const
 
template<class T >
void distribution2d (const std::vector< T > &inputX, const std::vector< T > &inputY, std::vector< std::vector< double > > &output, int nbin, T &minX, T &maxX, T &minY, T &maxY, double sigma=0, const std::string &filename="") const
 
template<class T >
void cumulative (const std::vector< T > &input, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end, std::vector< int > &output, int nbin, T &minimum, T &maximum) const
 
template<class T >
void percentiles (const std::vector< T > &input, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end, std::vector< T > &output, int nbin, T &minimum, T &maximum, const std::string &filename="") const
 
template<class T >
percentile (const std::vector< T > &input, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end, double percent, T minimum=0, T maximum=0) const
 
template<class T >
void signature (const std::vector< T > &input, double &k, double &alpha, double &beta, double e) const
 
void signature (double m1, double m2, double &k, double &alpha, double &beta, double e) const
 
template<class T >
void normalize (const std::vector< T > &input, std::vector< double > &output) const
 
template<class T >
void normalize_pct (std::vector< T > &input) const
 
template<class T >
double rmse (const std::vector< T > &x, const std::vector< T > &y) const
 
template<class T >
double nrmse (const std::vector< T > &x, const std::vector< T > &y) const
 
template<class T >
double cvrmse (const std::vector< T > &x, const std::vector< T > &y) const
 
template<class T >
double correlation (const std::vector< T > &x, const std::vector< T > &y, int delay=0) const
 
template<class T >
double gsl_covariance (const std::vector< T > &x, const std::vector< T > &y) const
 
template<class T >
double cross_correlation (const std::vector< T > &x, const std::vector< T > &y, int maxdelay, std::vector< T > &z) const
 
template<class T >
double linear_regression (const std::vector< T > &x, const std::vector< T > &y, double &c0, double &c1) const
 
template<class T >
double linear_regression_err (const std::vector< T > &x, const std::vector< T > &y, double &c0, double &c1) const
 
template<class T >
void interpolateNoData (const std::vector< double > &wavelengthIn, const std::vector< T > &input, const std::string &type, std::vector< T > &output, bool verbose=false) const
 
template<class T >
void interpolateUp (const std::vector< double > &wavelengthIn, const std::vector< T > &input, const std::vector< double > &wavelengthOut, const std::string &type, std::vector< T > &output, bool verbose=false) const
 
template<class T >
void interpolateUp (const std::vector< double > &wavelengthIn, const std::vector< std::vector< T > > &input, const std::vector< double > &wavelengthOut, const std::string &type, std::vector< std::vector< T > > &output, bool verbose=false) const
 
template<class T >
void interpolateUp (const std::vector< T > &input, std::vector< T > &output, int nbin) const
 
template<class T >
void nearUp (const std::vector< T > &input, std::vector< T > &output) const
 
template<class T >
void interpolateUp (double *input, int dim, std::vector< T > &output, int nbin)
 
template<class T >
void interpolateDown (const std::vector< T > &input, std::vector< T > &output, int nbin) const
 
template<class T >
void interpolateDown (double *input, int dim, std::vector< T > &output, int nbin)
 
template<class T >
std::vector< T >::const_iterator mymin (const std::vector< T > &v, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end) const
 
template<class T >
std::vector< T >::iterator mymin (const std::vector< T > &v, typename std::vector< T >::iterator begin, typename std::vector< T >::iterator end) const
 
template<class T >
std::vector< T >::const_iterator mymin (const std::vector< T > &v, typename std::vector< T >::const_iterator begin, typename std::vector< T >::const_iterator end, T minConstraint) const
 
template<class T >
std::vector< T >::iterator mymin (const std::vector< T > &v, typename std::vector< T >::iterator begin, typename std::vector< T >::iterator end, T minConstraint) const
 
template<class T >
mymin (const std::vector< T > &v) const
 
template<class T >
mymin (const std::vector< T > &v, T minConstraint) const
 
template<class T >
mymax (const std::vector< T > &v) const
 
template<class T >
mymax (const std::vector< T > &v, T maxConstraint) const
 
template<class T >
absmin (const std::vector< T > &v) const
 
template<class T >
absmax (const std::vector< T > &v) const
 

Static Public Member Functions

static void allocAcc (gsl_interp_accel *&acc)
 
static void getSpline (const std::string type, int size, gsl_spline *&spline)
 
static int initSpline (gsl_spline *spline, const double *x, const double *y, int size)
 
static double evalSpline (gsl_spline *spline, double x, gsl_interp_accel *acc)
 
static gsl_rng * getRandomGenerator (unsigned long int theSeed)
 

Detailed Description

Definition at line 43 of file StatFactory.h.

Member Enumeration Documentation

◆ DISTRIBUTION_TYPE

enum statfactory::StatFactory::DISTRIBUTION_TYPE

Definition at line 48 of file StatFactory.h.

48{uniform=1,gaussian=2};

◆ INTERPOLATION_TYPE

enum statfactory::StatFactory::INTERPOLATION_TYPE

Definition at line 46 of file StatFactory.h.

46{undefined=0,linear=1,polynomial=2,cspline=3,cspline_periodic=4,akima=5,akima_periodic=6};

Constructor & Destructor Documentation

◆ StatFactory()

statfactory::StatFactory::StatFactory ( void  )
inline

Definition at line 50 of file StatFactory.h.

50{};

◆ ~StatFactory()

virtual statfactory::StatFactory::~StatFactory ( void  )
inlinevirtual

Definition at line 51 of file StatFactory.h.

51{};

Member Function Documentation

◆ absmax() [1/2]

template<class T >
T statfactory::StatFactory::absmax ( const std::vector< T > &  v) const
inline

Definition at line 620 of file StatFactory.h.

621{
622 typename std::vector<T>::const_iterator tmpIt;
623 tmpIt=absmax(v, v.begin(), v.end());
624 if(tmpIt!=v.end())
625 return *tmpIt;
626 else if(m_noDataValues.size())
627 return m_noDataValues[0];
628 else{
629 std::string errorString="Error: no valid data found";
630 throw(errorString);
631 }
632}

◆ absmax() [2/2]

template<class T >
std::vector< T >::const_iterator statfactory::StatFactory::absmax ( const std::vector< T > &  v,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end 
) const
inline

Definition at line 598 of file StatFactory.h.

599{
600 bool isValid=false;
601 typename std::vector<T>::const_iterator tmpIt;
602 for (typename std::vector<T>::const_iterator it = begin; it!=end; ++it){
603 if(isNoData(*it))
604 continue;
605 if(isValid){
606 if(abs(*tmpIt)<abs(*it))
607 tmpIt=it;
608 }
609 else{
610 tmpIt=it;
611 isValid=true;
612 }
613 }
614 if(isValid)
615 return tmpIt;
616 else
617 return end;
618}

◆ absmin() [1/2]

template<class T >
T statfactory::StatFactory::absmin ( const std::vector< T > &  v) const
inline

Definition at line 584 of file StatFactory.h.

585{
586 typename std::vector<T>::const_iterator tmpIt;
587 tmpIt=absmin(v, v.begin(), v.end());
588 if(tmpIt!=v.end())
589 return *tmpIt;
590 else if(m_noDataValues.size())
591 return m_noDataValues[0];
592 else{
593 std::string errorString="Error: no valid data found";
594 throw(errorString);
595 }
596}

◆ absmin() [2/2]

template<class T >
std::vector< T >::const_iterator statfactory::StatFactory::absmin ( const std::vector< T > &  v,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end 
) const
inline

Definition at line 558 of file StatFactory.h.

559{
560 bool isValid=false;
561 if(v.empty()){
562 std::string errorString="Error: vector is empty";
563 throw(errorString);
564 }
565 typename std::vector<T>::const_iterator tmpIt;
566 for (typename std::vector<T>::const_iterator it = begin; it!=end; ++it){
567 if(isNoData(*it))
568 continue;
569 if(isValid){
570 if(abs(*tmpIt)>abs(*it))
571 tmpIt=it;
572 }
573 else{
574 tmpIt=it;
575 isValid=true;
576 }
577 }
578 if(isValid)
579 return tmpIt;
580 else
581 return end;
582}

◆ allocAcc()

static void statfactory::StatFactory::allocAcc ( gsl_interp_accel *&  acc)
inlinestatic

Definition at line 62 of file StatFactory.h.

62 {
63 acc = gsl_interp_accel_alloc ();
64 };

◆ cmoment()

template<class T >
double statfactory::StatFactory::cmoment ( const std::vector< T > &  v,
int  n 
) const

Definition at line 808 of file StatFactory.h.

809{
810 unsigned int validSize=0;
811 typename std::vector<T>::const_iterator it;
812 double m=0;
813 double m1=mean(v);
814 for(it = v.begin(); it!= v.end(); ++it){
815 if(isNoData(*it))
816 continue;
817 m+=pow((*it-m1),n);
818 ++validSize;
819 }
820 if(validSize)
821 return m/validSize;
822 else if(m_noDataValues.size())
823 return m_noDataValues[0];
824 else{
825 std::string errorString="Error: no valid data found";
826 throw(errorString);
827 }
828}

◆ correlation()

template<class T >
double statfactory::StatFactory::correlation ( const std::vector< T > &  x,
const std::vector< T > &  y,
int  delay = 0 
) const

Definition at line 1290 of file StatFactory.h.

1290 {
1291 double meanX=0;
1292 double meanY=0;
1293 double varX=0;
1294 double varY=0;
1295 double sXY=0;
1296 meanVar(x,meanX,varX);
1297 meanVar(y,meanY,varY);
1298 double denom = sqrt(varX*varY);
1299 bool isValid=false;
1300 if(denom){
1301 //Calculate the correlation series
1302 sXY = 0;
1303 for (int i=0;i<x.size();++i) {
1304 int j = i + delay;
1305 if (j < 0 || j >= y.size())
1306 continue;
1307 else if(isNoData(x[i])||isNoData(y[j]))
1308 continue;
1309 else{
1310 isValid=true;
1311 if(i<0){
1312 std::ostringstream s;
1313 s<<"Error: i must be positive";
1314 throw(s.str());
1315 }
1316 if(i>=x.size()){
1317 std::ostringstream s;
1318 s<<"Error: i must be smaller than x.size()";
1319 throw(s.str());
1320 }
1321 if(j<0){
1322 std::ostringstream s;
1323 s<<"Error: j must be positive";
1324 throw(s.str());
1325 }
1326 if(j>=y.size()){
1327 std::ostringstream s;
1328 s<<"Error: j must be smaller than y.size()";
1329 throw(s.str());
1330 }
1331 sXY += (x[i] - meanX) * (y[j] - meanY);
1332 }
1333 }
1334 if(isValid){
1335 double minSize=(x.size()<y.size())?x.size():y.size();
1336 return(sXY / denom / (minSize-1));
1337 }
1338 else if(m_noDataValues.size())
1339 return m_noDataValues[0];
1340 else{
1341 std::string errorString="Error: no valid data found";
1342 throw(errorString);
1343 }
1344 }
1345 else
1346 return 0;
1347}

◆ cross_correlation()

template<class T >
double statfactory::StatFactory::cross_correlation ( const std::vector< T > &  x,
const std::vector< T > &  y,
int  maxdelay,
std::vector< T > &  z 
) const

Definition at line 1350 of file StatFactory.h.

1350 {
1351 z.clear();
1352 double sumCorrelation=0;
1353 for (int delay=-maxdelay;delay<maxdelay;delay++) {
1354 z.push_back(correlation(x,y,delay));
1355 sumCorrelation+=z.back();
1356 }
1357 return sumCorrelation;
1358}

◆ cvrmse()

template<class T >
double statfactory::StatFactory::cvrmse ( const std::vector< T > &  x,
const std::vector< T > &  y 
) const

Definition at line 1255 of file StatFactory.h.

1255 {
1256 if(x.size()!=y.size()){
1257 std::ostringstream s;
1258 s<<"Error: x and y not equal in size";
1259 throw(s.str());
1260 }
1261 if(x.empty()){
1262 std::ostringstream s;
1263 s<<"Error: x is empty";
1264 throw(s.str());
1265 }
1266 std::vector<T> tmpX=x;
1267 eraseNoData(tmpX);
1268 std::vector<T> tmpY=y;
1269 eraseNoData(tmpY);
1270 double maxY=mymax(tmpY);
1271 double minY=mymin(tmpY);
1272 double rangeY=maxY-minY;
1273 double mse=0;
1274 for(int isample=0;isample<x.size();++isample){
1275 double e=x[isample]-y[isample];
1276 mse+=e*e/x.size();
1277 }
1278 return sqrt(mse)/mean(tmpY);
1279}

◆ distribution() [1/2]

template<class T >
void statfactory::StatFactory::distribution ( const std::vector< T > &  input,
std::vector< double > &  output,
int  nbin,
double  sigma = 0,
const std::string &  filename = "" 
) const
inline

Definition at line 174 of file StatFactory.h.

174{distribution(input,input.begin(),input.end(),output,nbin,0,0,sigma,filename);};

◆ distribution() [2/2]

template<class T >
void statfactory::StatFactory::distribution ( const std::vector< T > &  input,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end,
std::vector< double > &  output,
int  nbin,
T &  minimum,
T &  maximum,
double  sigma = 0,
const std::string &  filename = "" 
) const

Definition at line 889 of file StatFactory.h.

890{
891 double minValue=0;
892 double maxValue=0;
893 minmax(input,begin,end,minimum,maximum);
894 if(minimum<maximum&&minimum>minValue)
895 minValue=minimum;
896 if(minimum<maximum&&maximum<maxValue)
897 maxValue=maximum;
898
899 //todo: check...
900 minimum=minValue;
901 maximum=maxValue;
902
903 if(maximum<=minimum){
904 std::ostringstream s;
905 s<<"Error: could not calculate distribution (min>=max)";
906 throw(s.str());
907 }
908 if(!nbin){
909 std::string errorString="Error: nbin not defined";
910 throw(errorString);
911 }
912 if(!input.size()){
913 std::string errorString="Error: no valid data found";
914 throw(errorString);
915 }
916 if(output.size()!=nbin){
917 output.resize(nbin);
918 for(int i=0;i<nbin;output[i++]=0);
919 }
920 bool isValid=false;
921 typename std::vector<T>::const_iterator it;
922 for(it=begin;it!=end;++it){
923 if(*it<minimum)
924 continue;
925 if(*it>maximum)
926 continue;
927 if(isNoData(*it))
928 continue;
929 isValid=true;
930 if(sigma>0){
931 // minimum-=2*sigma;
932 // maximum+=2*sigma;
933 //create kde for Gaussian basis function
934 //todo: speed up by calculating first and last bin with non-zero contriubtion...
935 //todo: calculate real surface below pdf by using gsl_cdf_gaussian_P(x-mean+binsize,sigma)-gsl_cdf_gaussian_P(x-mean,sigma)
936 for(int ibin=0;ibin<nbin;++ibin){
937 double icenter=minimum+static_cast<double>(maximum-minimum)*(ibin+0.5)/nbin;
938 double thePdf=gsl_ran_gaussian_pdf(*it-icenter, sigma);
939 output[ibin]+=thePdf;
940 }
941 }
942 else{
943 int theBin=0;
944 if(*it==maximum)
945 theBin=nbin-1;
946 else if(*it>minimum && *it<maximum)
947 theBin=static_cast<int>(static_cast<double>((nbin-1)*(*it)-minimum)/(maximum-minimum));
948 ++output[theBin];
949 // if(*it==maximum)
950 // ++output[nbin-1];
951 // else if(*it>=minimum && *it<maximum)
952 // ++output[static_cast<int>(static_cast<double>((*it)-minimum)/(maximum-minimum)*nbin)];
953 }
954 }
955 if(!isValid){
956 std::string errorString="Error: no valid data found";
957 throw(errorString);
958 }
959 else if(!filename.empty()){
960 std::ofstream outputfile;
961 outputfile.open(filename.c_str());
962 if(!outputfile){
963 std::ostringstream s;
964 s<<"Error opening distribution file , " << filename;
965 throw(s.str());
966 }
967 for(int ibin=0;ibin<nbin;++ibin)
968 outputfile << minimum+static_cast<double>(maximum-minimum)*(ibin+0.5)/nbin << " " << static_cast<double>(output[ibin])/input.size() << std::endl;
969 outputfile.close();
970 }
971}

◆ distribution2d()

template<class T >
void statfactory::StatFactory::distribution2d ( const std::vector< T > &  inputX,
const std::vector< T > &  inputY,
std::vector< std::vector< double > > &  output,
int  nbin,
T &  minX,
T &  maxX,
T &  minY,
T &  maxY,
double  sigma = 0,
const std::string &  filename = "" 
) const

Definition at line 973 of file StatFactory.h.

974{
975 if(inputX.empty()){
976 std::ostringstream s;
977 s<<"Error: inputX is empty";
978 throw(s.str());
979 }
980 if(inputX.size()!=inputY.size()){
981 std::ostringstream s;
982 s<<"Error: inputX is empty";
983 throw(s.str());
984 }
985 unsigned long int npoint=inputX.size();
986 if(maxX<=minX)
987 minmax(inputX,inputX.begin(),inputX.end(),minX,maxX);
988 if(maxX<=minX){
989 std::ostringstream s;
990 s<<"Error: could not calculate distribution (minX>=maxX)";
991 throw(s.str());
992 }
993 if(maxY<=minY)
994 minmax(inputY,inputY.begin(),inputY.end(),minY,maxY);
995 if(maxY<=minY){
996 std::ostringstream s;
997 s<<"Error: could not calculate distribution (minY>=maxY)";
998 throw(s.str());
999 }
1000 if(nbin<=1){
1001 std::ostringstream s;
1002 s<<"Error: nbin must be larger than 1";
1003 throw(s.str());
1004 }
1005 output.resize(nbin);
1006 for(int i=0;i<nbin;++i){
1007 output[i].resize(nbin);
1008 for(int j=0;j<nbin;++j)
1009 output[i][j]=0;
1010 }
1011 int binX=0;
1012 int binY=0;
1013 for(unsigned long int ipoint=0;ipoint<npoint;++ipoint){
1014 if(inputX[ipoint]==maxX)
1015 binX=nbin-1;
1016 else
1017 binX=static_cast<int>(static_cast<double>(inputX[ipoint]-minX)/(maxX-minX)*nbin);
1018 if(inputY[ipoint]==maxY)
1019 binY=nbin-1;
1020 else
1021 binY=static_cast<int>(static_cast<double>(inputY[ipoint]-minY)/(maxY-minY)*nbin);
1022 if(binX<0){
1023 std::ostringstream s;
1024 s<<"Error: binX is smaller than 0";
1025 throw(s.str());
1026 }
1027 if(output.size()<=binX){
1028 std::ostringstream s;
1029 s<<"Error: output size must be larger than binX";
1030 throw(s.str());
1031 }
1032 if(binY<0){
1033 std::ostringstream s;
1034 s<<"Error: binY is smaller than 0";
1035 throw(s.str());
1036 }
1037 if(output.size()<=binY){
1038 std::ostringstream s;
1039 s<<"Error: output size must be larger than binY";
1040 throw(s.str());
1041 }
1042 if(sigma>0){
1043 // minX-=2*sigma;
1044 // maxX+=2*sigma;
1045 // minY-=2*sigma;
1046 // maxY+=2*sigma;
1047 //create kde for Gaussian basis function
1048 //todo: speed up by calculating first and last bin with non-zero contriubtion...
1049 for(int ibinX=0;ibinX<nbin;++ibinX){
1050 double centerX=minX+static_cast<double>(maxX-minX)*ibinX/nbin;
1051 double pdfX=gsl_ran_gaussian_pdf(inputX[ipoint]-centerX, sigma);
1052 for(int ibinY=0;ibinY<nbin;++ibinY){
1053 //calculate \integral_ibinX^(ibinX+1)
1054 double centerY=minY+static_cast<double>(maxY-minY)*ibinY/nbin;
1055 double pdfY=gsl_ran_gaussian_pdf(inputY[ipoint]-centerY, sigma);
1056 output[ibinX][binY]+=pdfX*pdfY;
1057 }
1058 }
1059 }
1060 else
1061 ++output[binX][binY];
1062 }
1063 if(!filename.empty()){
1064 std::ofstream outputfile;
1065 outputfile.open(filename.c_str());
1066 if(!outputfile){
1067 std::ostringstream s;
1068 s<<"Error opening distribution file , " << filename;
1069 throw(s.str());
1070 }
1071 for(int binX=0;binX<nbin;++binX){
1072 outputfile << std::endl;
1073 for(int binY=0;binY<nbin;++binY){
1074 double binValueX=0;
1075 if(nbin==maxX-minX+1)
1076 binValueX=minX+binX;
1077 else
1078 binValueX=minX+static_cast<double>(maxX-minX)*(binX+0.5)/nbin;
1079 double binValueY=0;
1080 if(nbin==maxY-minY+1)
1081 binValueY=minY+binY;
1082 else
1083 binValueY=minY+static_cast<double>(maxY-minY)*(binY+0.5)/nbin;
1084 double value=0;
1085 value=static_cast<double>(output[binX][binY])/npoint;
1086 outputfile << binValueX << " " << binValueY << " " << value << std::endl;
1087 /* double value=static_cast<double>(output[binX][binY])/npoint; */
1088 /* outputfile << (maxX-minX)*bin/(nbin-1)+minX << " " << (maxY-minY)*bin/(nbin-1)+minY << " " << value << std::endl; */
1089 }
1090 }
1091 outputfile.close();
1092 }
1093}

◆ eraseNoData()

template<class T >
void statfactory::StatFactory::eraseNoData ( std::vector< T > &  v) const
inline

Definition at line 721 of file StatFactory.h.

722{
723 if(m_noDataValues.size()){
724 typename std::vector<T>::iterator it=v.begin();
725 while(it!=v.end()){
726 if(isNoData(*it))
727 v.erase(it);
728 else
729 ++it;
730 }
731 }
732}

◆ evalSpline()

static double statfactory::StatFactory::evalSpline ( gsl_spline *  spline,
double  x,
gsl_interp_accel *  acc 
)
inlinestatic

Definition at line 95 of file StatFactory.h.

95 {
96 return gsl_spline_eval (spline, x, acc);
97 };

◆ getDistributionType()

DISTRIBUTION_TYPE statfactory::StatFactory::getDistributionType ( const std::string  distributionType)
inline

Definition at line 57 of file StatFactory.h.

57 {
58 std::map<std::string, DISTRIBUTION_TYPE> m_distMap;
59 initDist(m_distMap);
60 return m_distMap[distributionType];
61 };

◆ getInterpolationType()

INTERPOLATION_TYPE statfactory::StatFactory::getInterpolationType ( const std::string  interpolationType)
inline

Definition at line 52 of file StatFactory.h.

52 {
53 std::map<std::string, INTERPOLATION_TYPE> m_interpMap;
54 initMap(m_interpMap);
55 return m_interpMap[interpolationType];
56 };

◆ getNodataValues()

void statfactory::StatFactory::getNodataValues ( std::vector< double > &  nodatav) const
inline

Definition at line 108 of file StatFactory.h.

108{nodatav=m_noDataValues;};

◆ getRandomGenerator()

static gsl_rng * statfactory::StatFactory::getRandomGenerator ( unsigned long int  theSeed)
inlinestatic

Definition at line 99 of file StatFactory.h.

99 {
100 const gsl_rng_type * T;
101 gsl_rng * r;
102
103 // select random number generator
104 r = gsl_rng_alloc (gsl_rng_mt19937);
105 gsl_rng_set(r,theSeed);
106 return r;
107 }

◆ getRandomValue()

double statfactory::StatFactory::getRandomValue ( const gsl_rng *  r,
const std::string  type,
double  a = 0,
double  b = 1 
) const
inline

Definition at line 124 of file StatFactory.h.

124 {
125 std::map<std::string, DISTRIBUTION_TYPE> m_distMap;
126 initDist(m_distMap);
127 double randValue=0;
128 switch(m_distMap[type]){
129 case(uniform):
130 randValue = a+gsl_rng_uniform(r);
131 break;
132 case(gaussian):
133 default:
134 randValue = a+gsl_ran_gaussian(r, b);
135 break;
136 }
137 return randValue;
138 };

◆ getSpline()

static void statfactory::StatFactory::getSpline ( const std::string  type,
int  size,
gsl_spline *&  spline 
)
inlinestatic

Definition at line 66 of file StatFactory.h.

66 {
67 std::map<std::string, INTERPOLATION_TYPE> m_interpMap;
68 initMap(m_interpMap);
69 switch(m_interpMap[type]){
70 case(polynomial):
71 spline=gsl_spline_alloc(gsl_interp_polynomial,size);
72 break;
73 case(cspline):
74 spline=gsl_spline_alloc(gsl_interp_cspline,size);
75 break;
76 case(cspline_periodic):
77 spline=gsl_spline_alloc(gsl_interp_cspline_periodic,size);
78 break;
79 case(akima):
80 spline=gsl_spline_alloc(gsl_interp_akima,size);
81 break;
82 case(akima_periodic):
83 spline=gsl_spline_alloc(gsl_interp_akima_periodic,size);
84 break;
85 case(linear):
86 default:
87 spline=gsl_spline_alloc(gsl_interp_linear,size);
88 break;
89 }
90 assert(spline);
91 };

◆ gsl_covariance()

template<class T >
double statfactory::StatFactory::gsl_covariance ( const std::vector< T > &  x,
const std::vector< T > &  y 
) const

Definition at line 1285 of file StatFactory.h.

1285 {
1286 return(gsl_stats_covariance(&(x[0]),1,&(y[0]),1,x.size()));
1287 }

◆ initSpline()

static int statfactory::StatFactory::initSpline ( gsl_spline *  spline,
const double *  x,
const double *  y,
int  size 
)
inlinestatic

Definition at line 92 of file StatFactory.h.

92 {
93 return gsl_spline_init (spline, x, y, size);
94 };

◆ interpolateDown() [1/2]

template<class T >
void statfactory::StatFactory::interpolateDown ( const std::vector< T > &  input,
std::vector< T > &  output,
int  nbin 
) const

Definition at line 1608 of file StatFactory.h.

1609{
1610 if(input.empty()){
1611 std::ostringstream s;
1612 s<<"Error: input is empty";
1613 throw(s.str());
1614 }
1615 if(!nbin){
1616 std::ostringstream s;
1617 s<<"Error: nbin must be larger than 0";
1618 throw(s.str());
1619 }
1620 output.clear();
1621 int dim=input.size();
1622 int x=0;
1623 output.push_back(input[0]);
1624 for(int i=1;i<dim;++i){
1625 if(i%nbin)
1626 continue;
1627 else{
1628 x=(i-1)/nbin+1;
1629 output.push_back(input[i]);
1630 }
1631 }
1632}

◆ interpolateDown() [2/2]

template<class T >
void statfactory::StatFactory::interpolateDown ( double *  input,
int  dim,
std::vector< T > &  output,
int  nbin 
)

Definition at line 1635 of file StatFactory.h.

1636{
1637 if(!nbin){
1638 std::ostringstream s;
1639 s<<"Error: nbin must be larger than 0";
1640 throw(s.str());
1641 }
1642 output.clear();
1643 int x=0;
1644 output.push_back(input[0]);
1645 for(int i=1;i<dim;++i){
1646 if(i%nbin)
1647 continue;
1648 else{
1649 x=(i-1)/nbin+1;
1650 output.push_back(input[i]);
1651 }
1652 }
1653}

◆ interpolateNoData()

template<class T >
void statfactory::StatFactory::interpolateNoData ( const std::vector< double > &  wavelengthIn,
const std::vector< T > &  input,
const std::string &  type,
std::vector< T > &  output,
bool  verbose = false 
) const

Definition at line 1403 of file StatFactory.h.

1403 {
1404 if(wavelengthIn.empty()){
1405 std::ostringstream s;
1406 s<<"Error: wavelengthIn is empty";
1407 throw(s.str());
1408 }
1409 std::vector<double> wavelengthOut=wavelengthIn;
1410 std::vector<T> validIn=input;
1411 if(input.size()!=wavelengthIn.size()){
1412 std::ostringstream s;
1413 s<<"Error: x and y not equal in size";
1414 throw(s.str());
1415 }
1416 int nband=wavelengthIn.size();
1417 output.clear();
1418 //remove nodata from input and corresponding wavelengthIn
1419 if(m_noDataValues.size()){
1420 typename std::vector<T>::iterator itValue=validIn.begin();
1421 typename std::vector<T>::iterator itWavelength=wavelengthOut.begin();
1422 while(itValue!=validIn.end()&&itWavelength!=wavelengthOut.end()){
1423 if(isNoData(*itValue)){
1424 validIn.erase(itValue);
1425 wavelengthOut.erase(itWavelength);
1426 }
1427 else{
1428 ++itValue;
1429 ++itWavelength;
1430 }
1431 }
1432 if(validIn.size()>1){
1433 try{
1434 interpolateUp(wavelengthOut, validIn, wavelengthIn, type, output, verbose);
1435 }
1436 catch(...){
1437 output=input;
1438 }
1439 }
1440 else//we can not interpolate if no valid data
1441 output=input;
1442 }
1443 else//no nodata values to interpolate
1444 output=input;
1445}

◆ interpolateUp() [1/3]

template<class T >
void statfactory::StatFactory::interpolateUp ( const std::vector< double > &  wavelengthIn,
const std::vector< T > &  input,
const std::vector< double > &  wavelengthOut,
const std::string &  type,
std::vector< T > &  output,
bool  verbose = false 
) const

Definition at line 1447 of file StatFactory.h.

1447 {
1448 if(wavelengthIn.empty()){
1449 std::ostringstream s;
1450 s<<"Error: wavelengthIn is empty";
1451 throw(s.str());
1452 }
1453 if(input.size()!=wavelengthIn.size()){
1454 std::ostringstream s;
1455 s<<"Error: input and wavelengthIn not equal in size";
1456 throw(s.str());
1457 }
1458 if(wavelengthOut.empty()){
1459 std::ostringstream s;
1460 s<<"Error: wavelengthOut is empty";
1461 throw(s.str());
1462 }
1463 int nband=wavelengthIn.size();
1464 output.clear();
1465 gsl_interp_accel *acc;
1466 allocAcc(acc);
1467 gsl_spline *spline;
1468 getSpline(type,nband,spline);
1469 assert(spline);
1470 assert(&(wavelengthIn[0]));
1471 assert(&(input[0]));
1472 int status=initSpline(spline,&(wavelengthIn[0]),&(input[0]),nband);
1473 if(status){
1474 std::string errorString="Could not initialize spline";
1475 throw(errorString);
1476 }
1477 for(int index=0;index<wavelengthOut.size();++index){
1478 if(wavelengthOut[index]<*wavelengthIn.begin()){
1479 output.push_back(*(input.begin()));
1480 continue;
1481 }
1482 else if(wavelengthOut[index]>wavelengthIn.back()){
1483 output.push_back(input.back());
1484 continue;
1485 }
1486 double dout=evalSpline(spline,wavelengthOut[index],acc);
1487 output.push_back(dout);
1488 }
1489 gsl_spline_free(spline);
1490 gsl_interp_accel_free(acc);
1491}

◆ interpolateUp() [2/3]

template<class T >
void statfactory::StatFactory::interpolateUp ( const std::vector< T > &  input,
std::vector< T > &  output,
int  nbin 
) const

Definition at line 1525 of file StatFactory.h.

1526{
1527 if(input.empty()){
1528 std::ostringstream s;
1529 s<<"Error: input is empty";
1530 throw(s.str());
1531 }
1532 if(!nbin){
1533 std::ostringstream s;
1534 s<<"Error: nbin must be larger than 0";
1535 throw(s.str());
1536 }
1537 output.clear();
1538 int dim=input.size();
1539 for(int i=0;i<dim;++i){
1540 double deltaX=0;
1541 double left=input[i];
1542 if(i<dim-1){
1543 double right=(i<dim-1)? input[i+1]:input[i];
1544 deltaX=(right-left)/static_cast<double>(nbin);
1545 for(int x=0;x<nbin;++x){
1546 output.push_back(left+x*deltaX);
1547 }
1548 }
1549 else
1550 output.push_back(input.back());
1551 }
1552}

◆ interpolateUp() [3/3]

template<class T >
void statfactory::StatFactory::interpolateUp ( double *  input,
int  dim,
std::vector< T > &  output,
int  nbin 
)

Definition at line 1584 of file StatFactory.h.

1585{
1586 if(!nbin){
1587 std::ostringstream s;
1588 s<<"Error: nbin must be larger than 0";
1589 throw(s.str());
1590 }
1591 output.clear();
1592 for(int i=0;i<dim;++i){
1593 double deltaX=0;
1594 double left=input[i];
1595 if(i<dim-1){
1596 double right=(i<dim-1)? input[i+1]:input[i];
1597 deltaX=(right-left)/static_cast<double>(nbin);
1598 for(int x=0;x<nbin;++x){
1599 output.push_back(left+x*deltaX);
1600 }
1601 }
1602 else
1603 output.push_back(input[dim-1]);
1604 }
1605}

◆ isNoData()

bool statfactory::StatFactory::isNoData ( double  value) const
inline

Definition at line 109 of file StatFactory.h.

109 {
110 if(m_noDataValues.empty())
111 return false;
112 else
113 return find(m_noDataValues.begin(),m_noDataValues.end(),value)!=m_noDataValues.end();
114 };

◆ kurtosis()

template<class T >
double statfactory::StatFactory::kurtosis ( const std::vector< T > &  v) const

Definition at line 836 of file StatFactory.h.

837{
838 //todo: what if nodata value?
839 double m2=cmoment(v,2);
840 double m4=cmoment(v,4);
841 return m4/m2/m2-3.0;
842}

◆ linear_regression()

template<class T >
double statfactory::StatFactory::linear_regression ( const std::vector< T > &  x,
const std::vector< T > &  y,
double &  c0,
double &  c1 
) const

Definition at line 1361 of file StatFactory.h.

1361 {
1362 if(x.size()!=y.size()){
1363 std::ostringstream s;
1364 s<<"Error: x and y not equal in size";
1365 throw(s.str());
1366 }
1367 if(x.empty()){
1368 std::ostringstream s;
1369 s<<"Error: x is empty";
1370 throw(s.str());
1371 }
1372 double cov00;
1373 double cov01;
1374 double cov11;
1375 double sumsq;
1376 gsl_fit_linear(&(x[0]),1,&(y[0]),1,x.size(),&c0,&c1,&cov00,&cov01,&cov11,&sumsq);
1377 return (1-sumsq/var(y)/(y.size()-1));
1378}

◆ linear_regression_err()

template<class T >
double statfactory::StatFactory::linear_regression_err ( const std::vector< T > &  x,
const std::vector< T > &  y,
double &  c0,
double &  c1 
) const

Definition at line 1381 of file StatFactory.h.

1381 {
1382 if(x.size()!=y.size()){
1383 std::ostringstream s;
1384 s<<"Error: x and y not equal in size";
1385 throw(s.str());
1386 }
1387 if(x.empty()){
1388 std::ostringstream s;
1389 s<<"Error: x is empty";
1390 throw(s.str());
1391 }
1392 double cov00;
1393 double cov01;
1394 double cov11;
1395 double sumsq;
1396 gsl_fit_linear(&(x[0]),1,&(y[0]),1,x.size(),&c0,&c1,&cov00,&cov01,&cov11,&sumsq);
1397 return sqrt((sumsq)/(y.size()));
1398}

◆ mean()

template<class T >
double statfactory::StatFactory::mean ( const std::vector< T > &  v) const
inline

Definition at line 700 of file StatFactory.h.

701{
702 typename std::vector<T>::const_iterator it;
703 T tmpSum=0;
704 unsigned int validSize=0;
705 for (it = v.begin(); it!= v.end(); ++it){
706 if(isNoData(*it))
707 continue;
708 ++validSize;
709 tmpSum+=*it;
710 }
711 if(validSize)
712 return static_cast<double>(tmpSum)/validSize;
713 else if(m_noDataValues.size())
714 return m_noDataValues[0];
715 else{
716 std::string errorString="Error: no valid data found";
717 throw(errorString);
718 }
719}

◆ meanVar()

template<class T >
void statfactory::StatFactory::meanVar ( const std::vector< T > &  v,
double &  m1,
double &  v1 
) const

Definition at line 844 of file StatFactory.h.

845{
846 typename std::vector<T>::const_iterator it;
847 unsigned int validSize=0;
848 m1=0;
849 v1=0;
850 double m2=0;
851 for (it = v.begin(); it!= v.end(); ++it){
852 if(isNoData(*it))
853 continue;
854 m1+=*it;
855 m2+=(*it)*(*it);
856 ++validSize;
857 }
858 if(validSize){
859 m2/=validSize;
860 m1/=validSize;
861 v1=m2-m1*m1;
862 }
863 else if(m_noDataValues.size()){
864 m1=m_noDataValues[0];
865 v1=m_noDataValues[0];
866 }
867 else{
868 std::string errorString="Error: no valid data found";
869 throw(errorString);
870 }
871}

◆ median()

template<class T >
T statfactory::StatFactory::median ( const std::vector< T > &  v) const

Definition at line 740 of file StatFactory.h.

741{
742 std::vector<T> tmpV=v;
743 eraseNoData(tmpV);
744 if(tmpV.size()){
745 sort(tmpV.begin(),tmpV.end());
746 if(tmpV.size()%2)
747 return tmpV[tmpV.size()/2];
748 else
749 return 0.5*(tmpV[tmpV.size()/2-1]+tmpV[tmpV.size()/2]);
750 }
751 else if(m_noDataValues.size())
752 return m_noDataValues[0];
753 else{
754 std::string errorString="Error: no valid data found";
755 throw(errorString);
756 }
757}

◆ minmax()

template<class T >
void statfactory::StatFactory::minmax ( const std::vector< T > &  v,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end,
T &  theMin,
T &  theMax 
) const
inline

Definition at line 634 of file StatFactory.h.

635{
636 bool isConstraint=(theMax>theMin);
637 double minConstraint=theMin;
638 double maxConstraint=theMax;
639 bool isValid=false;
640 for (typename std::vector<T>::const_iterator it = begin; it!=end; ++it){
641 if(isNoData(*it))
642 continue;
643 if(isValid){
644 if(isConstraint){
645 if(*it<minConstraint)
646 continue;
647 if(*it>maxConstraint)
648 continue;
649 }
650 if(*it<theMin)
651 theMin=*it;
652 if(*it>theMax)
653 theMax=*it;
654 }
655 else{
656 if(isConstraint){
657 if(*it<minConstraint)
658 continue;
659 if(*it>maxConstraint)
660 continue;
661 }
662 theMin=*it;
663 theMax=*it;
664 isValid=true;
665 }
666 }
667 if(!isValid){
668 if(m_noDataValues.size()){
669 theMin=m_noDataValues[0];
670 theMax=m_noDataValues[0];
671 }
672 else{
673 std::string errorString="Error: no valid data found";
674 throw(errorString);
675 }
676 }
677}

◆ moment()

template<class T >
double statfactory::StatFactory::moment ( const std::vector< T > &  v,
int  n 
) const

Definition at line 785 of file StatFactory.h.

786{
787 unsigned int validSize=0;
788 typename std::vector<T>::const_iterator it;
789 double m=0;
790// double m1=mean(v);
791 for(it = v.begin(); it!= v.end(); ++it){
792 if(isNoData(*it))
793 continue;
794 m+=pow((*it),n);
795 ++validSize;
796 }
797 if(validSize)
798 return m/validSize;
799 else if(m_noDataValues.size())
800 return m_noDataValues[0];
801 else{
802 std::string errorString="Error: no valid data found";
803 throw(errorString);
804 }
805}

◆ mymax() [1/6]

template<class T >
T statfactory::StatFactory::mymax ( const std::vector< T > &  v) const
inline

Definition at line 500 of file StatFactory.h.

501{
502 bool isValid=false;
503 if(v.empty()){
504 std::string errorString="Error: vector is empty";
505 throw(errorString);
506 }
507 T maxValue;
508 for (typename std::vector<T>::const_iterator it = v.begin(); it!=v.end(); ++it){
509 if(isNoData(*it))
510 continue;
511 if(isValid){
512 if(maxValue<*it)
513 maxValue=*it;
514 }
515 else{
516 maxValue=*it;
517 isValid=true;
518 }
519 }
520 if(isValid)
521 return maxValue;
522 else if(m_noDataValues.size())
523 return m_noDataValues[0];
524 else{
525 std::string errorString="Error: no valid data found";
526 throw(errorString);
527 }
528}

◆ mymax() [2/6]

template<class T >
T statfactory::StatFactory::mymax ( const std::vector< T > &  v,
maxConstraint 
) const
inline

Definition at line 530 of file StatFactory.h.

531{
532 bool isValid=false;
533 T maxValue=maxConstraint;
534 for (typename std::vector<T>::const_iterator it = v.begin(); it!=v.end(); ++it){
535 if(isNoData(*it))
536 continue;
537 if(isValid){
538 if((*it<=maxConstraint)&&(*it>maxValue))
539 maxValue=*it;
540 }
541 else{
542 if(*it>maxValue)
543 continue;
544 maxValue=*it;
545 isValid=true;
546 }
547 }
548 if(isValid)
549 return maxValue;
550 else if(m_noDataValues.size())
551 return m_noDataValues[0];
552 else{
553 std::string errorString="Error: no valid data found";
554 throw(errorString);
555 }
556}

◆ mymax() [3/6]

template<class T >
std::vector< T >::const_iterator statfactory::StatFactory::mymax ( const std::vector< T > &  v,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end 
) const
inline

Definition at line 338 of file StatFactory.h.

339{
340 bool isValid=false;
341 typename std::vector<T>::const_iterator tmpIt;
342 for (typename std::vector<T>::iterator it = begin; it!=end; ++it){
343 if(isNoData(*it))
344 continue;
345 if(isValid){
346 if(*tmpIt<*it)
347 tmpIt=it;
348 }
349 else{
350 tmpIt=it;
351 isValid=true;
352 }
353 }
354 if(isValid)
355 return tmpIt;
356 else if(m_noDataValues.size())
357 return m_noDataValues[0];
358 else{
359 std::string errorString="Error: no valid data found";
360 throw(errorString);
361 }
362}

◆ mymax() [4/6]

template<class T >
std::vector< T >::const_iterator statfactory::StatFactory::mymax ( const std::vector< T > &  v,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end,
maxConstraint 
) const
inline

Definition at line 386 of file StatFactory.h.

387{
388 bool isValid=false;
389 typename std::vector<T>::const_iterator tmpIt;
390 T maxValue=maxConstraint;
391 for (typename std::vector<T>::const_iterator it = begin; it!=end; ++it){
392 if(isNoData(*it))
393 continue;
394 if(isValid){
395 if((maxConstraint>=*it)&&(*it>maxValue)){
396 tmpIt=it;
397 maxValue=*it;
398 }
399 }
400 else{
401 if(*it>maxConstraint)
402 continue;
403 tmpIt=it;
404 maxValue=*it;
405 isValid=true;
406 }
407 }
408 if(isValid)
409 return tmpIt;
410 else
411 return end;
412}

◆ mymax() [5/6]

template<class T >
std::vector< T >::iterator statfactory::StatFactory::mymax ( const std::vector< T > &  v,
typename std::vector< T >::iterator  begin,
typename std::vector< T >::iterator  end 
) const
inline

Definition at line 364 of file StatFactory.h.

365{
366 bool isValid=false;
367 typename std::vector<T>::iterator tmpIt;
368 for (typename std::vector<T>::iterator it = begin; it!=end; ++it){
369 if(isNoData(*it))
370 continue;
371 if(isValid){
372 if(*tmpIt<*it)
373 tmpIt=it;
374 }
375 else{
376 tmpIt=it;
377 isValid=true;
378 }
379 }
380 if(isValid)
381 return tmpIt;
382 else
383 return end;
384}

◆ mymax() [6/6]

template<class T >
std::vector< T >::iterator statfactory::StatFactory::mymax ( const std::vector< T > &  v,
typename std::vector< T >::iterator  begin,
typename std::vector< T >::iterator  end,
maxConstraint 
) const
inline

Definition at line 414 of file StatFactory.h.

415{
416 bool isValid=false;
417 typename std::vector<T>::iterator tmpIt=v.end();
418 T maxValue=maxConstraint;
419 for (typename std::vector<T>::const_iterator it = begin; it!=end; ++it){
420 if(isNoData(*it))
421 continue;
422 if(isValid){
423 if((maxConstraint>=*it)&&(*it>maxValue)){
424 tmpIt=it;
425 maxValue=*it;
426 }
427 }
428 else{
429 if(*it>maxValue)
430 continue;
431 tmpIt=it;
432 maxValue=*it;
433 isValid=true;
434 }
435 }
436 if(isValid)
437 return tmpIt;
438 else
439 return end;
440}

◆ mymin() [1/6]

template<class T >
T statfactory::StatFactory::mymin ( const std::vector< T > &  v) const
inline

Definition at line 442 of file StatFactory.h.

443{
444 bool isValid=false;
445 if(v.empty()){
446 std::string errorString="Error: vector is empty";
447 throw(errorString);
448 }
449 T minValue;
450 for (typename std::vector<T>::const_iterator it = v.begin(); it!=v.end(); ++it){
451 if(isNoData(*it))
452 continue;
453 if(isValid){
454 if(minValue>*it)
455 minValue=*it;
456 }
457 else{
458 minValue=*it;
459 isValid=true;
460 }
461 }
462 if(isValid)
463 return minValue;
464 else if(m_noDataValues.size())
465 return m_noDataValues[0];
466 else{
467 std::string errorString="Error: no valid data found";
468 throw(errorString);
469 }
470}

◆ mymin() [2/6]

template<class T >
T statfactory::StatFactory::mymin ( const std::vector< T > &  v,
minConstraint 
) const
inline

Definition at line 472 of file StatFactory.h.

473{
474 bool isValid=false;
475 T minValue=minConstraint;
476 for (typename std::vector<T>::const_iterator it = v.begin(); it!=v.end(); ++it){
477 if(isNoData(*it))
478 continue;
479 if(isValid){
480 if((minConstraint<=*it)&&(*it<minValue))
481 minValue=*it;
482 }
483 else{
484 if(*it<minValue)
485 continue;
486 minValue=*it;
487 isValid=true;
488 }
489 }
490 if(isValid)
491 return minValue;
492 else if(m_noDataValues.size())
493 return m_noDataValues[0];
494 else{
495 std::string errorString="Error: no valid data found";
496 throw(errorString);
497 }
498}

◆ mymin() [3/6]

template<class T >
std::vector< T >::const_iterator statfactory::StatFactory::mymin ( const std::vector< T > &  v,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end 
) const
inline

Definition at line 222 of file StatFactory.h.

223{
224 bool isValid=false;
225 typename std::vector<T>::const_iterator tmpIt;
226 for(typename std::vector<T>::const_iterator it = begin; it!=end; ++it){
227 if(!isNoData(*it)){
228 if(isValid){
229 if(*tmpIt>*it)
230 tmpIt=it;
231 }
232 else{
233 tmpIt=it;
234 isValid=true;
235 }
236 }
237 }
238 if(isValid)
239 return tmpIt;
240 else if(m_noDataValues.size())
241 return m_noDataValues[0];
242 else{
243 std::string errorString="Error: no valid data found";
244 throw(errorString);
245 }
246}

◆ mymin() [4/6]

template<class T >
std::vector< T >::const_iterator statfactory::StatFactory::mymin ( const std::vector< T > &  v,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end,
minConstraint 
) const
inline

Definition at line 274 of file StatFactory.h.

275{
276 bool isValid=false;
277 typename std::vector<T>::const_iterator tmpIt;
278 T minValue=minConstraint;
279 for (typename std::vector<T>::const_iterator it = begin; it!=end; ++it){
280 if(isNoData(*it))
281 continue;
282 if(isValid){
283 if((minConstraint<=*it)&&(*it<minValue)){
284 tmpIt=it;
285 minValue=*it;
286 }
287 }
288 else{
289 if(*it<minValue)
290 continue;
291 tmpIt=it;
292 minValue=*it;
293 isValid=true;
294 }
295 }
296 if(isValid)
297 return tmpIt;
298 else if(m_noDataValues.size())
299 return m_noDataValues[0];
300 else{
301 std::string errorString="Error: no valid data found";
302 throw(errorString);
303 }
304}

◆ mymin() [5/6]

template<class T >
std::vector< T >::iterator statfactory::StatFactory::mymin ( const std::vector< T > &  v,
typename std::vector< T >::iterator  begin,
typename std::vector< T >::iterator  end 
) const
inline

Definition at line 248 of file StatFactory.h.

249{
250 bool isValid=false;
251 typename std::vector<T>::iterator tmpIt;
252 for (typename std::vector<T>::const_iterator it = begin; it!=end; ++it){
253 if(!isNoData(*it)){
254 if(isValid){
255 if(*tmpIt>*it)
256 tmpIt=it;
257 }
258 else{
259 tmpIt=it;
260 isValid=true;
261 }
262 }
263 }
264 if(isValid)
265 return tmpIt;
266 else if(m_noDataValues.size())
267 return m_noDataValues[0];
268 else{
269 std::string errorString="Error: no valid data found";
270 throw(errorString);
271 }
272}

◆ mymin() [6/6]

template<class T >
std::vector< T >::iterator statfactory::StatFactory::mymin ( const std::vector< T > &  v,
typename std::vector< T >::iterator  begin,
typename std::vector< T >::iterator  end,
minConstraint 
) const
inline

Definition at line 306 of file StatFactory.h.

307{
308 bool isValid=false;
309 typename std::vector<T>::iterator tmpIt;
310 T minValue=minConstraint;
311 for (typename std::vector<T>::const_iterator it = begin; it!=end; ++it){
312 if(isNoData(*it))
313 continue;
314 if(isValid){
315 if((minConstraint<=*it)&&(*it<minValue)){
316 tmpIt=it;
317 minValue=*it;
318 }
319 }
320 else{
321 if(*it<minConstraint)
322 continue;
323 tmpIt=it;
324 minValue=*it;
325 isValid=true;
326 }
327 }
328 if(isValid)
329 return tmpIt;
330 else if(m_noDataValues.size())
331 return m_noDataValues[0];
332 else{
333 std::string errorString="Error: no valid data found";
334 throw(errorString);
335 }
336}

◆ nearUp()

template<class T >
void statfactory::StatFactory::nearUp ( const std::vector< T > &  input,
std::vector< T > &  output 
) const

Definition at line 1555 of file StatFactory.h.

1556{
1557 if(input.empty()){
1558 std::ostringstream s;
1559 s<<"Error: input is empty";
1560 throw(s.str());
1561 }
1562 if(output.size()<input.size()){
1563 std::ostringstream s;
1564 s<<"Error: output size is smaller than input size";
1565 throw(s.str());
1566 }
1567 int dimInput=input.size();
1568 int dimOutput=output.size();
1569
1570 for(int iin=0;iin<dimInput;++iin){
1571 for(int iout=0;iout<dimOutput/dimInput;++iout){
1572 int indexOutput=iin*dimOutput/dimInput+iout;
1573 if(indexOutput>=output.size()){
1574 std::ostringstream s;
1575 s<<"Error: indexOutput must be smaller than output.size()";
1576 throw(s.str());
1577 }
1578 output[indexOutput]=input[iin];
1579 }
1580 }
1581}

◆ normalize()

template<class T >
void statfactory::StatFactory::normalize ( const std::vector< T > &  input,
std::vector< double > &  output 
) const

Definition at line 1185 of file StatFactory.h.

1185 {
1186 double total=sum(input);
1187 if(total){
1188 output.resize(input.size());
1189 for(int index=0;index<input.size();++index)
1190 output[index]=input[index]/total;
1191 }
1192 else
1193 output=input;
1194}

◆ normalize_pct()

template<class T >
void statfactory::StatFactory::normalize_pct ( std::vector< T > &  input) const

Definition at line 1197 of file StatFactory.h.

1197 {
1198 double total=sum(input);
1199 if(total){
1200 typename std::vector<T>::iterator it;
1201 for(it=input.begin();it!=input.end();++it)
1202 *it=100.0*(*it)/total;
1203 }
1204}

◆ nrmse()

template<class T >
double statfactory::StatFactory::nrmse ( const std::vector< T > &  x,
const std::vector< T > &  y 
) const

Definition at line 1228 of file StatFactory.h.

1228 {
1229 if(x.size()!=y.size()){
1230 std::ostringstream s;
1231 s<<"Error: x and y not equal in size";
1232 throw(s.str());
1233 }
1234 if(x.empty()){
1235 std::ostringstream s;
1236 s<<"Error: x is empty";
1237 throw(s.str());
1238 }
1239 std::vector<T> tmpX=x;
1240 eraseNoData(tmpX);
1241 std::vector<T> tmpY=y;
1242 eraseNoData(tmpY);
1243 double maxY=mymax(y);
1244 double minY=mymin(y);
1245 double rangeY=maxY-minY;
1246 double mse=0;
1247 for(int isample=0;isample<x.size();++isample){
1248 double e=x[isample]-y[isample];
1249 mse+=e*e/x.size();
1250 }
1251 return sqrt(mse)/rangeY;
1252}

◆ nvalid()

template<class T >
unsigned int statfactory::StatFactory::nvalid ( const std::vector< T > &  v) const

Definition at line 734 of file StatFactory.h.

734 {
735 std::vector<T> tmpV=v;
736 eraseNoData(tmpV);
737 return(tmpV.size());
738 }

◆ percentile()

template<class T >
T statfactory::StatFactory::percentile ( const std::vector< T > &  input,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end,
double  percent,
minimum = 0,
maximum = 0 
) const

Definition at line 1154 of file StatFactory.h.

1155{
1156 if(input.empty()){
1157 std::ostringstream s;
1158 s<<"Error: input is empty";
1159 throw(s.str());
1160 }
1161 std::vector<double> inputSort;
1162 inputSort.assign(begin,end);
1163 typename std::vector<double>::iterator vit=inputSort.begin();
1164 while(vit!=inputSort.end()){
1165 if(maximum>minimum){
1166 if(*vit<minimum||*vit>maximum)
1167 inputSort.erase(vit);
1168 }
1169 else
1170 ++vit;
1171 }
1172 eraseNoData(inputSort);
1173 std::sort(inputSort.begin(),inputSort.end());
1174 return gsl_stats_quantile_from_sorted_data(&(inputSort[0]),1,inputSort.size(),percent/100.0);
1175}

◆ percentiles()

template<class T >
void statfactory::StatFactory::percentiles ( const std::vector< T > &  input,
typename std::vector< T >::const_iterator  begin,
typename std::vector< T >::const_iterator  end,
std::vector< T > &  output,
int  nbin,
T &  minimum,
T &  maximum,
const std::string &  filename = "" 
) const

Definition at line 1096 of file StatFactory.h.

1097{
1098 if(maximum<=minimum)
1099 minmax(input,begin,end,minimum,maximum);
1100 if(maximum<=minimum){
1101 std::ostringstream s;
1102 s<<"Error: maximum must be at least minimum";
1103 throw(s.str());
1104 }
1105 if(nbin<=1){
1106 std::ostringstream s;
1107 s<<"Error: nbin must be larger than 1";
1108 throw(s.str());
1109 }
1110 if(input.empty()){
1111 std::ostringstream s;
1112 s<<"Error: input is empty";
1113 throw(s.str());
1114 }
1115 output.resize(nbin);
1116 std::vector<T> inputSort;
1117 inputSort.assign(begin,end);
1118 typename std::vector<T>::iterator vit=inputSort.begin();
1119 while(vit!=inputSort.end()){
1120 if(*vit<minimum||*vit>maximum)
1121 inputSort.erase(vit);
1122 else
1123 ++vit;
1124 }
1125 eraseNoData(inputSort);
1126 std::sort(inputSort.begin(),inputSort.end());
1127 vit=inputSort.begin();
1128 std::vector<T> inputBin;
1129 for(int ibin=0;ibin<nbin;++ibin){
1130 inputBin.clear();
1131 while(inputBin.size()<inputSort.size()/nbin&&vit!=inputSort.end()){
1132 inputBin.push_back(*vit);
1133 ++vit;
1134 }
1135 if(inputBin.size()){
1136 output[ibin]=mymax(inputBin);
1137 }
1138 }
1139 if(!filename.empty()){
1140 std::ofstream outputfile;
1141 outputfile.open(filename.c_str());
1142 if(!outputfile){
1143 std::ostringstream s;
1144 s<<"error opening distribution file , " << filename;
1145 throw(s.str());
1146 }
1147 for(int ibin=0;ibin<nbin;++ibin)
1148 outputfile << ibin*100.0/nbin << " " << static_cast<double>(output[ibin])/input.size() << std::endl;
1149 outputfile.close();
1150 }
1151}

◆ pushNodDataValue()

unsigned int statfactory::StatFactory::pushNodDataValue ( double  noDataValue)
inline

Definition at line 115 of file StatFactory.h.

115 {
116 if(find(m_noDataValues.begin(),m_noDataValues.end(),noDataValue)==m_noDataValues.end())
117 m_noDataValues.push_back(noDataValue);
118 return m_noDataValues.size();
119 };

◆ rmse()

template<class T >
double statfactory::StatFactory::rmse ( const std::vector< T > &  x,
const std::vector< T > &  y 
) const

Definition at line 1206 of file StatFactory.h.

1206 {
1207 if(x.size()!=y.size()){
1208 std::ostringstream s;
1209 s<<"Error: x and y not equal in size";
1210 throw(s.str());
1211 }
1212 if(x.empty()){
1213 std::ostringstream s;
1214 s<<"Error: x is empty";
1215 throw(s.str());
1216 }
1217 double mse=0;
1218 for(int isample=0;isample<x.size();++isample){
1219 if(isNoData(x[isample])||isNoData(y[isample]))
1220 continue;
1221 double e=x[isample]-y[isample];
1222 mse+=e*e/x.size();
1223 }
1224 return sqrt(mse);
1225}

◆ scale2byte()

template<class T1 , class T2 >
void statfactory::StatFactory::scale2byte ( const std::vector< T1 > &  input,
std::vector< T2 > &  output,
unsigned char  lbound = 0,
unsigned char  ubound = 255 
) const

Definition at line 873 of file StatFactory.h.

874{
875 output.resize(input.size());
876 T1 minimum=mymin(input);
877 T1 maximum=mymax(input);
878 if(minimum>=maximum){
879 std::string errorString="Error: no valid data found";
880 throw(errorString);
881 }
882 double scale=(ubound-lbound)/(maximum-minimum);
883 //todo: what if nodata value?
884 for (int i=0;i<input.size();++i){
885 output[i]=scale*(input[i]-(minimum))+lbound;
886 }
887}

◆ setNoDataValues()

unsigned int statfactory::StatFactory::setNoDataValues ( std::vector< double >  vnodata)
inline

Definition at line 120 of file StatFactory.h.

120 {
121 m_noDataValues=vnodata;
122 return m_noDataValues.size();
123 };

◆ signature()

template<class T >
void statfactory::StatFactory::signature ( const std::vector< T > &  input,
double &  k,
double &  alpha,
double &  beta,
double  e 
) const

Definition at line 1177 of file StatFactory.h.

1178{
1179 double m1=moment(input,1);
1180 double m2=moment(input,2);
1181 signature(m1,m2,k,alpha,beta,e);
1182}

◆ skewness()

template<class T >
double statfactory::StatFactory::skewness ( const std::vector< T > &  v) const

Definition at line 830 of file StatFactory.h.

831{
832 //todo: what if nodata value?
833 return cmoment(v,3)/pow(var(v),1.5);
834}

◆ sum()

template<class T >
T statfactory::StatFactory::sum ( const std::vector< T > &  v) const
inline

Definition at line 679 of file StatFactory.h.

680{
681 bool isValid=false;
682 typename std::vector<T>::const_iterator it;
683 T tmpSum=0;
684 for (it = v.begin(); it!= v.end(); ++it){
685 if(isNoData(*it))
686 continue;
687 isValid=true;
688 tmpSum+=*it;
689 }
690 if(isValid)
691 return tmpSum;
692 else if(m_noDataValues.size())
693 return m_noDataValues[0];
694 else{
695 std::string errorString="Error: no valid data found";
696 throw(errorString);
697 }
698}

◆ var()

template<class T >
double statfactory::StatFactory::var ( const std::vector< T > &  v) const

Definition at line 759 of file StatFactory.h.

760{
761 typename std::vector<T>::const_iterator it;
762 unsigned int validSize=0;
763 double m1=0;
764 double m2=0;
765 for (it = v.begin(); it!= v.end(); ++it){
766 if(isNoData(*it))
767 continue;
768 m1+=*it;
769 m2+=(*it)*(*it);
770 ++validSize;
771 }
772 if(validSize){
773 m2/=validSize;
774 m1/=validSize;
775 return m2-m1*m1;
776 }
777 else if(m_noDataValues.size())
778 return m_noDataValues[0];
779 else{
780 std::string errorString="Error: no valid data found";
781 throw(errorString);
782 }
783}

The documentation for this class was generated from the following file: