| pktools 2.6.7 Processing Kernel for geospatial data

Vector2d.h
1/**********************************************************************
2 Vector2d.h: 2-dimensional vector class (inherits from stl vector class)
3 Copyright (C) 2008-2012 Pieter Kempeneers
4
5 This file is part of pktools
6
7 pktools is free software: you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10(at your option) any later version.
11
12 pktools is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with pktools. If not, see <http://www.gnu.org/licenses/>.
19***********************************************************************/
20#ifndef _VECTOR2D_H_
21#define _VECTOR2D_H_
22
23#include <vector>
24#include <list>
25#include <algorithm>
26#include <numeric>
27#include <gsl/gsl_matrix.h>
28#include "IndexValue.h"
29#include "algorithms/StatFactory.h"
30
31 template<class T> class Vector2d: public std::vector<std::vector <T> >
32{
33 public:
34 Vector2d();
35 Vector2d(const Vector2d<T>& v1);//copy constructor
36 ~Vector2d();
37 Vector2d(int nrow);
38 Vector2d(int nrow, int ncol);
39 Vector2d(int nrow, int ncol, const T& value);
40 Vector2d(const gsl_matrix* gsl_m);
41 void resize(int nrow)
42 {
43 std::vector< std::vector<T> >::resize(nrow);
44 };
45 void resize(int nrow, int ncol);
46 int nRows() const {return this->size();};
47 int nCols() const {if(this->size()) return this->begin()->size(); else return 0;};
48 void selectCol(int col, std::vector<T> &output) const;
49 void selectCol(int col, T* output) const;
50 std::vector<T> selectCol(int col);
51 void selectCols(const std::list<int> &cols, Vector2d<T> &output) const;
52 void setMask(const Vector2d<T> &mask, T msknodata, T nodata=0);
53 void transpose(Vector2d<T> &output) const{
54 output.resize(nCols(),nRows());
55 for(int irow=0;irow<nRows();++irow){
56 for(int icol=0;icol<nCols();++icol){
57 output[icol][irow]=(*this)[irow][icol];
58 }
59 }
60 };
61 void selectCols(const std::list<int> &cols);
62 void sort(Vector2d<T>& output);
63 void scale(const std::vector<double> &scaleVector, const std::vector<double> &offsetVector, Vector2d<T>& scaledOutput);
64 void scale(const T lbound, const T ubound, std::vector<double> &scaleVector, std::vector<double> &offsetVector, Vector2d<T>& scaledOutput);
65 Vector2d<T> operator=(const Vector2d<T>& v1);
66 Vector2d<T> operator+=(const Vector2d<T>& v1);
67// std::ostream& operator<<(std::ostream& os, const Vector2d<T>& v);
68// template<class T> std::ostream& operator<<(std::ostream& os, const Vector2d<T>& v);
69 template<class T1> friend std::ostream& operator<<(std::ostream & os, const Vector2d<T1>& v);
70 Vector2d<T> sum(const Vector2d<T>& v1, const Vector2d<T>& v2) const;
71 T mymax(int& x, int& y, double maxValue) const;
72
73 T sum() const;
74};
75
76 template<class T> Vector2d<T>::Vector2d()
77 : std::vector< std::vector<T> >()
78{
79}
80
81 template<class T> Vector2d<T>::~Vector2d()
82{
83}
84
85//copy constructor
86 template<class T> Vector2d<T>::Vector2d(const Vector2d<T>& v1){
87 this->resize(v1.size());
88 for(int irow=0;irow<v1.size();++irow)
89 this->at(irow)=v1[irow];
90}
91
92 template<class T> Vector2d<T> Vector2d<T>::operator=(const Vector2d<T>& v1){
93 //check for assignment to self (of the form v=v)
94 if(this==&v1)
95 return *this;
96 else{
97 this->resize(v1.size());
98 for(int irow=0;irow<v1.size();++irow)
99 this->at(irow)=v1[irow];
100 return *this;
101 }
102}
103
104 template<class T> Vector2d<T> Vector2d<T>::operator+=(const Vector2d<T>& v1){
105 assert(v1.nRows()==nRows());
106 assert(v1.nCols()==nCols());
107 for(int irow=0;irow<nRows();++irow)
108 for(int icol=0;icol<nCols();++icol)
109 (*this)[irow][icol]+=v1[irow][icol];
110 return *this;
111}
112
113 template<class T> Vector2d<T>::Vector2d(int nrow)
114 : std::vector< std::vector<T> >(nrow)
115{
116}
117
118 template<class T> Vector2d<T>::Vector2d(int nrow, int ncol)
119// : std::vector< std::vector<T> >(nrow)
120{
121 this->resize(nrow);
122 for(int irow=0;irow<nrow;++irow){
123 (this->operator[](irow)).resize(ncol);
124// (*this)[irow].resize(ncol);
125 }
126}
127
128 template<class T> Vector2d<T>::Vector2d(int nrow, int ncol, const T& value)
129{
130 this->resize(nrow);
131 for(int irow=0;irow<nrow;++irow){
132 (this->operator[](irow)).resize(ncol);
133 for(int icol=0;icol<ncol;++icol)
134 (this->operator[](irow))[icol]=value;
135 }
136}
137
138 template<class T> Vector2d<T>::Vector2d(const gsl_matrix* gsl_m)
139{
140 this->resize(gsl_m->size1);
141 for(int irow=0;irow<this->size();++irow){
142 (this->operator[](irow)).resize(gsl_m->size2);
143 for(int icol=0;icol<this->operator[](irow).size();++icol)
144 (this->operator[](irow))[icol]=gsl_matrix_get(gsl_m,irow,icol);
145 }
146}
147
148
149 template<class T> void Vector2d<T>::resize(int nrow, int ncol)
150{
151 this->std::vector< std::vector<T> >::resize(nrow);
152 for(int irow=0;irow<nrow;++irow){
153 (this->operator[](irow)).resize(ncol);
154 }
155}
156
157 template<class T> void Vector2d<T>::selectCols(const std::list<int> &cols, Vector2d<T> &output) const
158{
159 output.resize(this->size());
160 std::list<int>::const_iterator it;
161 for(int irow=0;irow<this->size();++irow){
162 output[irow].resize(cols.size());
163 it=cols.begin();
164 for(int icol=0;icol<cols.size();++icol)
165 output[irow][icol]=(*this)[irow][*(it++)];
166 }
167}
168
169 template<class T> void Vector2d<T>::selectCol(int col, std::vector<T> &output) const
170{
171 assert(col>=0);
172 assert(col<(*this)[0].size());
173 output.resize(this->size());
174 for(int irow=0;irow<this->size();++irow){
175 output[irow]=(*this)[irow][col];
176 }
177}
178
179 template<class T> std::vector<T> Vector2d<T>::selectCol(int col)
180{
181 assert(col>=0);
182 assert(col<(*this)[0].size());
183 std::vector<T> output(this->size());
184 for(int irow=0;irow<this->size();++irow)
185 output[irow]=(*this)[irow][col];
186 return(output);
187}
188
189 template<class T> void Vector2d<T>::selectCol(int col, T* output) const
190{
191 assert(col>=0);
192 assert(col<(*this)[0].size());
193 for(int irow=0;irow<this->size();++irow){
194 output[irow]=(*this)[irow][col];
195 }
196}
197
198 template<class T> void Vector2d<T>::selectCols(const std::list<int> &cols)
199{
200 for(int irow=0;irow<this->size();++irow)
201 for(int icol=((*this)[irow]).size()-1;icol>=0;--icol)
202 if(find(cols.begin(),cols.end(),icol)==cols.end())
203 (*this)[irow].erase(((*this)[irow]).begin()+icol);
204}
205
206 template<class T> void Vector2d<T>::setMask(const Vector2d<T> &mask, T msknodata, T nodata)
207{
208 assert(mask.nRows()==nRows());
209 assert(mask.nCols()==nCols());
210 for(int irow=0;irow<this->size();++irow)
211 for(int icol=0;icol<((*this)[irow]).size()-1;++icol)
212 if(mask[irow][icol]==msknodata)
213 (*this)[irow][icol]=nodata;
214}
215
216 template<class T1> std::ostream& operator<<(std::ostream& os, const Vector2d<T1>& v)
217{
218 for(int irow=0;irow<v.size();++irow){
219 for(int icol=0;icol<v[irow].size();++icol){
220 os << v[irow][icol] << "\t";
221 }
222 os << std::endl;
223 }
224 return os;
225 // os << theOption.getLongName() << ": ";
226 // for(int index=0;index<theOption.size();++index)
227 // os << type2string<T>(theOption[index]) << " ";
228 // os << std::endl;
229 // return os;
230}
231
232 template<class T> void Vector2d<T>::sort(Vector2d<T>& output)
233{
234 //sort according to first sample (ex. wavelength)
235 int nsample=this->size();//including first sample (ex. wavelength)
236 int nband=(*this)[0].size();
237 std::vector<IndexValue> sortW(nband);
238 for(int ilevel=0;ilevel<nband;++ilevel){
239 IndexValue pv;
240 pv.position=ilevel;
241 pv.value=(*this)[0][ilevel];
242 sortW[ilevel]=pv;
243 }
244 std::sort(sortW.begin(),sortW.end(),Increase_IndexValue());
245 output.resize(nsample);
246 for(int isample=0;isample<nsample;++isample){
247 output[isample].resize(nband);
248 for(int iband=0;iband<nband;++iband)
249 output[isample][iband]=(*this)[isample][sortW[iband].position];
250 }
251}
252
253 template<class T> void Vector2d<T>::scale(const std::vector<double> &scaleVector,const std::vector<double> &offsetVector, Vector2d<T>& scaledOutput)
254{
255 int nsample=this->size();//including first sample (ex. wavelength)
256 int nband=(*this)[0].size();
257 assert(scaleVector.size()==nband);
258 assert(offsetVector.size()==nband);
259 std::vector<T> pixel(nband);
260 scaledOutput.resize(nsample,nband);
261 for(int isample=0;isample<nsample;++isample)
262 for(int iband=0;iband<nband;++iband)
263 scaledOutput[isample][iband]=((*this)[isample][iband])*scaleVector[iband]+offsetVector[iband];
264}
265
266 template<class T> void Vector2d<T>::scale(const T lbound, const T ubound, std::vector<double> &scaleVector, std::vector<double> &offsetVector, Vector2d<T>& scaledOutput)
267{
268 //scale to lbound and ubound
269 int nsample=this->size();//including first sample (ex. wavelength)
270 int nband=(*this)[0].size();
271 scaleVector.resize(nband);
272 offsetVector.resize(nband);
273 std::vector<T> pixel(nsample);
274 T theMin;
275 T theMax;
276 statfactory::StatFactory stat;
277 scaledOutput.resize(nsample,nband);
278 for(int iband=0;iband<nband;++iband){
279 pixel=selectCol(iband);
280 stat.minmax(pixel, pixel.begin(), pixel.end(), theMin, theMax);
281 scaleVector[iband]=static_cast<double>(ubound-lbound)/(theMax-theMin);
282 offsetVector[iband]=static_cast<double>(-theMin*scaleVector[iband])-lbound;
283 for(int isample=0;isample<pixel.size();++isample)
284 scaledOutput[isample][iband]=((*this)[isample][iband])*scaleVector[iband]+offsetVector[iband];
285 }
286}
287
288 template<class T> Vector2d<T> Vector2d<T>::sum(const Vector2d<T>& v1, const Vector2d<T>& v2) const{
289 Vector2d<T> vsum(v1.size());
290 assert(v1.size()==v2.size());
291 for(int irow=0;irow<v1.size();++irow){
292 assert(v1[irow].size()==v2[irow].size());
293 vsum[irow].resize(v1[irow].size());
294 for(int icol=0;icol<v1.size();++icol)
295 vsum[irow][icol]=v1[irow][icol]+v2[irow][icol];
296 }
297 return vsum;
298}
299
300 template<class T> T Vector2d<T>::sum() const{
301 double theSum=0;
302 for(int irow=0;irow<this->size();++irow){
303 for(int icol=0;icol<this->operator[](irow).size();++icol)
304 theSum+=(this->operator[](irow))[icol];
305 }
306 return theSum;
307}
308
309 template<class T> T Vector2d<T>::mymax(int& x, int& y, double maxValue) const{
310 //todo: what if this->operator[](0)[0] >=maxValue?
311 // double theMax=(this->operator[](0))[0];
312 double theMax=0;
313 for(int irow=0;irow<this->size();++irow){
314 for(int icol=0;icol<(this->operator[](irow)).size();++icol){
315 double currentValue=(this->operator[](irow))[icol];
316 if(currentValue<maxValue&&currentValue>theMax){
317 assert(theMax<maxValue);
318 y=irow;
319 x=icol;
320 theMax=currentValue;
321 }
322 }
323 }
324 assert(theMax<maxValue);
325 return theMax;
326}
327
328#endif /* _VECTOR2D_H_ */
Increase_IndexValue
Definition: IndexValue.h:39
Vector2d
Definition: Vector2d.h:32
statfactory::StatFactory
Definition: StatFactory.h:43
IndexValue
Definition: IndexValue.h:23