This is an automated email from the git hooks/post-receive script. sebastic-guest pushed a commit to branch upstream-master in repository pktools.
commit 8e25ca92a1e4f2751965103729316bd7ae94ab9a Author: Pieter Kempeneers <kempe...@gmail.com> Date: Thu Apr 25 17:38:35 2013 +0200 removed dependency of gslwrap: replaced gslwrap with armadillo in SensorModel.h and pksensormodel, replaced gslwrap with stl::vector and Vector2d --- src/algorithms/Filter.h | 34 ++---- src/algorithms/Makefile.am | 2 +- src/apps/Makefile.am | 5 +- src/apps/pksensormodel.cc | 182 +++++++++++++-------------- src/apps/pksensormodel.h | 25 ++-- src/models/SensorModel.h | 298 ++++++++++++++++++++++----------------------- 6 files changed, 265 insertions(+), 281 deletions(-) diff --git a/src/algorithms/Filter.h b/src/algorithms/Filter.h index 9d2eb9d..3da2601 100644 --- a/src/algorithms/Filter.h +++ b/src/algorithms/Filter.h @@ -22,9 +22,6 @@ along with pktools. If not, see <http://www.gnu.org/licenses/>. #include <vector> #include <iostream> -#include <gslwrap/vector_double.h> -#include <gslwrap/matrix_double.h> -#include <gslwrap/matrix_vector_operators.h> #include "StatFactory.h" #include "imageclasses/ImgReaderGdal.h" #include "imageclasses/ImgWriterGdal.h" @@ -120,9 +117,6 @@ private: gsl_spline *spline; stat.getSpline(interpolationType,nband,spline); stat.initSpline(spline,&(srf[0][0]),&(srf[1][0]),nband); - // gsl_interp_accel *acc=gsl_interp_accel_alloc(); - // gsl_spline *spline=gsl_spline_alloc(gsl_interp_linear,nband); - // gsl_spline_init(spline,&(srf[0][0]),&(srf[1][0]),nband); if(verbose) std::cout << "calculating norm of srf" << std::endl << std::flush; double norm=0; @@ -201,9 +195,6 @@ private: gsl_spline *spline; stat.getSpline(interpolationType,nband,spline); stat.initSpline(spline,&(srf[0][0]),&(srf[1][0]),nband); - // gsl_interp_accel *acc=gsl_interp_accel_alloc(); - // gsl_spline *spline=gsl_spline_alloc(gsl_interp_linear,nband); - // gsl_spline_init(spline,&(srf[0][0]),&(srf[1][0]),nband); if(verbose) std::cout << "calculating norm of srf" << std::endl << std::flush; double norm=0; @@ -292,22 +283,23 @@ template<class T> void Filter::applyFwhm(const vector<double> &wavelengthIn, con int nbandOut=wavelengthOut.size(); output.resize(nbandOut); - gsl::matrix tf(nbandIn,nbandOut); + Vector2d<double> tf(nbandIn,nbandOut); for(int indexOut=0;indexOut<nbandOut;++indexOut){ double norm=0; for(int indexIn=0;indexIn<nbandIn;++indexIn){ - tf(indexIn,indexOut)= + // tf(indexIn,indexOut)= + tf[indexIn][indexOut]= exp((wavelengthOut[indexOut]-wavelength_fine[indexIn]) *(wavelength_fine[indexIn]-wavelengthOut[indexOut]) /2.0/stddev[indexOut] /stddev[indexOut]); - tf(indexIn,indexOut)/=sqrt(2.0*M_PI); - tf(indexIn,indexOut)/=stddev[indexOut]; - norm+=tf(indexIn,indexOut); + tf[indexIn][indexOut]/=sqrt(2.0*M_PI); + tf[indexIn][indexOut]/=stddev[indexOut]; + norm+=tf[indexIn][indexOut]; } output[indexOut]=0; for(int indexIn=0;indexIn<nbandIn;++indexIn) - output[indexOut]+=input_fine[indexIn]*tf(indexIn,indexOut)/norm; + output[indexOut]+=input_fine[indexIn]*tf[indexIn][indexOut]/norm; } } @@ -335,19 +327,19 @@ template<class T> void Filter::applyFwhm(const vector<double> &wavelengthIn, con int nbandOut=wavelengthOut.size(); output.resize(nbandOut,(input[0].size()+down-1)/down); - gsl::matrix tf(nbandIn,nbandOut); + Vector2d<double> tf(nbandIn,nbandOut); vector<double> norm(nbandOut); for(int indexOut=0;indexOut<nbandOut;++indexOut){ norm[indexOut]=0; for(int indexIn=0;indexIn<nbandIn;++indexIn){ - tf(indexIn,indexOut)= + tf[indexIn][indexOut]= exp((wavelengthOut[indexOut]-wavelength_fine[indexIn]) *(wavelength_fine[indexIn]-wavelengthOut[indexOut]) /2.0/stddev[indexOut] /stddev[indexOut]); - tf(indexIn,indexOut)/=sqrt(2.0*M_PI); - tf(indexIn,indexOut)/=stddev[indexOut]; - norm[indexOut]+=tf(indexIn,indexOut); + tf[indexIn][indexOut]/=sqrt(2.0*M_PI); + tf[indexIn][indexOut]/=stddev[indexOut]; + norm[indexOut]+=tf[indexIn][indexOut]; } } @@ -362,7 +354,7 @@ template<class T> void Filter::applyFwhm(const vector<double> &wavelengthIn, con stat.interpolateUp(wavelengthIn,inputValues,wavelength_fine,interpolationType,input_fine,verbose); output[indexOut][(isample+down-1)/down]=0; for(int indexIn=0;indexIn<nbandIn;++indexIn){ - output[indexOut][(isample+down-1)/down]+=input_fine[indexIn]*tf(indexIn,indexOut)/norm[indexOut]; + output[indexOut][(isample+down-1)/down]+=input_fine[indexIn]*tf[indexIn][indexOut]/norm[indexOut]; } } } diff --git a/src/algorithms/Makefile.am b/src/algorithms/Makefile.am index d9b35a9..b3f9028 100644 --- a/src/algorithms/Makefile.am +++ b/src/algorithms/Makefile.am @@ -37,4 +37,4 @@ libalgorithms_la_SOURCES = $(libalgorithms_la_HEADERS) Egcs.cc Filter2d.cc Filte # list of sources for the binaries #pktestStat_SOURCES = pktestStat.cc -#pktestStat_LDADD = $(GSL_LIBS) -lgslwrap $(top_builddir)/src/algorithms/libalgorithms.la $(top_builddir)/src/imageclasses/libimageClasses.la -lgslwrap +#pktestStat_LDADD = $(GSL_LIBS) $(top_builddir)/src/algorithms/libalgorithms.la $(top_builddir)/src/imageclasses/libimageClasses.la diff --git a/src/apps/Makefile.am b/src/apps/Makefile.am index b64c749..f80400e 100644 --- a/src/apps/Makefile.am +++ b/src/apps/Makefile.am @@ -29,8 +29,9 @@ pklas2img_LDADD = -llas $(AM_LDFLAGS) endif if USE_NLOPT bin_PROGRAMS += pksensormodel pkopt_svm +pksensormodel_LDFLAGS = -O1 pksensormodel_SOURCES = $(top_srcdir)/src/models/SensorModel.h pksensormodel.h pksensormodel.cc -pksensormodel_LDADD = $(GSL_LIBS) $(AM_LDFLAGS) -lnlopt -lm -lgslwrap +pksensormodel_LDADD = $(GSL_LIBS) $(AM_LDFLAGS) -lnlopt -lm -larmadillo #-lgslwrap pkopt_svm_SOURCES = $(top_srcdir)/src/algorithms/OptFactory.h pkclassify_nn.h pkopt_svm.cc pkopt_svm_LDADD = $(GSL_LIBS) $(AM_LDFLAGS) -lnlopt endif @@ -53,7 +54,7 @@ pkegcs_LDADD = -lgdal $(AM_LDFLAGS) -lgdal pkextract_SOURCES = pkextract.cc pkfillnodata_SOURCES = pkfillnodata.cc pkfilter_SOURCES = pkfilter.cc -pkfilter_LDADD = $(GSL_LIBS) $(AM_LDFLAGS) -lgsl -lgslwrap -lgdal +pkfilter_LDADD = $(GSL_LIBS) $(AM_LDFLAGS) -lgsl -lgdal #-lgslwrap pkdsm2shadow_SOURCES = pkdsm2shadow.cc pkmosaic_SOURCES = pkmosaic.cc pkndvi_SOURCES = pkndvi.cc diff --git a/src/apps/pksensormodel.cc b/src/apps/pksensormodel.cc index a63a589..24594c7 100644 --- a/src/apps/pksensormodel.cc +++ b/src/apps/pksensormodel.cc @@ -32,10 +32,10 @@ double objFunction(const std::vector<double> &x, std::vector<double> &grad, void assert(grad.empty()); double error=0; DataModel *dm=reinterpret_cast<DataModel*> (my_func_data); - gsl::vector bc_att(3); - bc_att[0]=x[0]; - bc_att[1]=x[1]; - bc_att[2]=x[2]; + arma::vec bc_att(3); + bc_att(0)=x[0]; + bc_att(1)=x[1]; + bc_att(2)=x[2]; dm->setBoresightAtt(bc_att); for(unsigned int ipoint=0;ipoint<dm->getSize();++ipoint){ double e=dm->getDistGeo(ipoint); @@ -190,8 +190,8 @@ int main(int argc, char *argv[]) theModel.setPPy(ppy_opt[0]); theModel.setPolynome(polynome_opt); theModel.setDatum(datum_opt[0]); - // gsl::vector bc_pos(3); - gsl::vector bc_att(3); + // arma::vec bc_pos(3); + arma::vec bc_att(3); // while(bcpos_opt.size()<3) // bcpos_opt.push_back(bcpos_opt[0]); while(bcatt_opt.size()<3) @@ -199,16 +199,16 @@ int main(int argc, char *argv[]) if(bcrad_opt[0]){ // bc_pos[0]=theModel.rad2deg(bcpos_opt[0]); // bc_pos[1]=theModel.rad2deg(bcpos_opt[1]); - bc_att[0]=theModel.rad2deg(bcatt_opt[0]); - bc_att[1]=theModel.rad2deg(bcatt_opt[1]); - bc_att[2]=theModel.rad2deg(bcatt_opt[2]); + bc_att(0)=theModel.rad2deg(bcatt_opt[0]); + bc_att(1)=theModel.rad2deg(bcatt_opt[1]); + bc_att(2)=theModel.rad2deg(bcatt_opt[2]); } else{ - // bc_pos[0]=bcpos_opt[0]; - // bc_pos[1]=bcpos_opt[1]; - bc_att[0]=bcatt_opt[0]; - bc_att[1]=bcatt_opt[1]; - bc_att[2]=bcatt_opt[2]; + // bc_pos(0)=bcpos_opt(0); + // bc_pos(1)=bcpos_opt(1); + bc_att(0)=bcatt_opt[0]; + bc_att(1)=bcatt_opt[1]; + bc_att(2)=bcatt_opt[2]; } // bc_pos[2]=bcpos_opt[2]; // theModel.setBoresightPos(bc_pos); @@ -246,9 +246,9 @@ int main(int argc, char *argv[]) istringstream csvstream(csvRecord); string item; int icol=0;//first column is id - gsl::vector thePosGCP(3); - gsl::vector thePosPlatform(3); - gsl::vector theAttPlatform(3); + arma::vec thePosGCP(3); + arma::vec thePosPlatform(3); + arma::vec theAttPlatform(3); while(getline(csvstream,item,fs_opt[0])){//read a column if(verbose_opt[0]>1) std::cout << item << " "; @@ -268,59 +268,59 @@ int main(int argc, char *argv[]) case(3)://X dvalue=atof(item.c_str()); if(gcprad_opt[0]) - thePosGCP[0]=theModel.rad2deg(dvalue); + thePosGCP(0)=theModel.rad2deg(dvalue); else - thePosGCP[0]=dvalue; + thePosGCP(0)=dvalue; break; case(4)://Y dvalue=atof(item.c_str()); if(gcprad_opt[0]) - thePosGCP[1]=theModel.rad2deg(dvalue); + thePosGCP(1)=theModel.rad2deg(dvalue); else - thePosGCP[1]=dvalue; + thePosGCP(1)=dvalue; break; case(5)://Z dvalue=atof(item.c_str()); - thePosGCP[2]=dvalue; + thePosGCP(2)=dvalue; break; case(6)://XL dvalue=atof(item.c_str()); if(pplrad_opt[0]) - thePosPlatform[0]=theModel.rad2deg(dvalue); + thePosPlatform(0)=theModel.rad2deg(dvalue); else - thePosPlatform[0]=dvalue; + thePosPlatform(0)=dvalue; break; case(7)://YL dvalue=atof(item.c_str()); if(pplrad_opt[0]) - thePosPlatform[1]=theModel.rad2deg(dvalue); + thePosPlatform(1)=theModel.rad2deg(dvalue); else - thePosPlatform[1]=dvalue; + thePosPlatform(1)=dvalue; break; case(8)://ZL dvalue=atof(item.c_str()); - thePosPlatform[2]=dvalue; + thePosPlatform(2)=dvalue; break; case(9)://Roll dvalue=atof(item.c_str()); if(aplrad_opt[0]) - theAttPlatform[0]=theModel.rad2deg(dvalue); + theAttPlatform(0)=theModel.rad2deg(dvalue); else - theAttPlatform[0]=dvalue; + theAttPlatform(0)=dvalue; break; case(10)://Pitch dvalue=atof(item.c_str()); if(aplrad_opt[0]) - theAttPlatform[1]=theModel.rad2deg(dvalue); + theAttPlatform(1)=theModel.rad2deg(dvalue); else - theAttPlatform[1]=dvalue; + theAttPlatform(1)=dvalue; break; case(11)://Yaw dvalue=atof(item.c_str()); if(aplrad_opt[0]) - theAttPlatform[2]=theModel.rad2deg(dvalue); + theAttPlatform(2)=theModel.rad2deg(dvalue); else - theAttPlatform[2]=dvalue; + theAttPlatform(2)=dvalue; break; case(12)://track break; @@ -332,10 +332,10 @@ int main(int argc, char *argv[]) } if(s_srs_opt[0]!=4326){ if(verbose_opt[0]>1) - std::cout << "transforming: " << thePosGCP[0] << ", " << thePosGCP[1] << ", " << thePosGCP[2] << std::endl; - poCT->Transform(1,&(thePosGCP[0]),&(thePosGCP[1]),&(thePosGCP[2])); + std::cout << "transforming: " << thePosGCP(0) << ", " << thePosGCP(1) << ", " << thePosGCP(2) << std::endl; + poCT->Transform(1,&(thePosGCP(0)),&(thePosGCP(1)),&(thePosGCP(2))); if(verbose_opt[0]>1) - std::cout << "into: " << thePosGCP[0] << ", " << thePosGCP[1] << ", " << thePosGCP[2] << std::endl; + std::cout << "into: " << thePosGCP(0) << ", " << thePosGCP(1) << ", " << thePosGCP(2) << std::endl; // poCT->Transform(1,&(thePosPlatform[0]),&(thePosPlatform[1]),&(thePosPlatform[2])); } theDataModel.pushPosGCP(thePosGCP); @@ -354,9 +354,9 @@ int main(int argc, char *argv[]) istringstream lineStream(spaceRecord); string item; int icol=0; - gsl::vector thePosGCP(3); - gsl::vector thePosPlatform(3); - gsl::vector theAttPlatform(3); + arma::vec thePosGCP(3); + arma::vec thePosPlatform(3); + arma::vec theAttPlatform(3); while(lineStream >> item){ if(verbose_opt[0]>1) std::cout << item << " "; @@ -377,59 +377,59 @@ int main(int argc, char *argv[]) case(3)://X itemStream >> dvalue; if(gcprad_opt[0]) - thePosGCP[0]=theModel.rad2deg(dvalue); + thePosGCP(0)=theModel.rad2deg(dvalue); else - thePosGCP[0]=dvalue; + thePosGCP(0)=dvalue; break; case(4)://Y itemStream >> dvalue; if(gcprad_opt[0]) - thePosGCP[1]=theModel.rad2deg(dvalue); + thePosGCP(1)=theModel.rad2deg(dvalue); else - thePosGCP[1]=dvalue; + thePosGCP(1)=dvalue; break; case(5)://Z itemStream >> dvalue; - thePosGCP[2]=dvalue; + thePosGCP(2)=dvalue; break; case(6)://XL itemStream >> dvalue; if(pplrad_opt[0]) - thePosPlatform[0]=theModel.rad2deg(dvalue); + thePosPlatform(0)=theModel.rad2deg(dvalue); else - thePosPlatform[0]=dvalue; + thePosPlatform(0)=dvalue; break; case(7)://YL itemStream >> dvalue; if(pplrad_opt[0]) - thePosPlatform[1]=theModel.rad2deg(dvalue); + thePosPlatform(1)=theModel.rad2deg(dvalue); else - thePosPlatform[1]=dvalue; + thePosPlatform(1)=dvalue; break; case(8)://ZL itemStream >> dvalue; - thePosPlatform[2]=dvalue; + thePosPlatform(2)=dvalue; break; case(9)://Roll itemStream >> dvalue; if(aplrad_opt[0]) - theAttPlatform[0]=theModel.rad2deg(dvalue); + theAttPlatform(0)=theModel.rad2deg(dvalue); else - theAttPlatform[0]=dvalue; + theAttPlatform(0)=dvalue; break; case(10)://Pitch itemStream >> dvalue; if(aplrad_opt[0]) - theAttPlatform[1]=theModel.rad2deg(dvalue); + theAttPlatform(1)=theModel.rad2deg(dvalue); else - theAttPlatform[1]=dvalue; + theAttPlatform(1)=dvalue; break; case(11)://Yaw itemStream >> dvalue; if(aplrad_opt[0]) - theAttPlatform[2]=theModel.rad2deg(dvalue); + theAttPlatform(2)=theModel.rad2deg(dvalue); else - theAttPlatform[2]=dvalue; + theAttPlatform(2)=dvalue; break; case(12)://track break; @@ -441,10 +441,10 @@ int main(int argc, char *argv[]) } if(s_srs_opt[0]!=4326){ if(verbose_opt[0]>1) - std::cout << "transforming: " << thePosGCP[0] << ", " << thePosGCP[1] << ", " << thePosGCP[2] << std::endl; - poCT->Transform(1,&(thePosGCP[0]),&(thePosGCP[1]),&(thePosGCP[2])); + std::cout << "transforming: " << thePosGCP(0) << ", " << thePosGCP(1) << ", " << thePosGCP(2) << std::endl; + poCT->Transform(1,&(thePosGCP(0)),&(thePosGCP(1)),&(thePosGCP(2))); if(verbose_opt[0]>1) - std::cout << "into: " << thePosGCP[0] << ", " << thePosGCP[1] << ", " << thePosGCP[2] << std::endl; + std::cout << "into: " << thePosGCP(0) << ", " << thePosGCP(1) << ", " << thePosGCP(2) << std::endl; // poCT->Transform(1,&(thePosPlatform[0]),&(thePosPlatform[1]),&(thePosPlatform[2])); } theDataModel.pushPosGCP(thePosGCP); @@ -465,8 +465,8 @@ int main(int argc, char *argv[]) if(verbose_opt[0]>1) std::cout << "ipoint: " << ipoint << std::endl; vID.push_back(ipoint); - gsl::vector thePosPlatform(3); - gsl::vector theAttPlatform(3); + arma::vec thePosPlatform(3); + arma::vec theAttPlatform(3); assert(row_opt.size()>ipoint); assert(col_opt.size()>ipoint); theDataModel.pushRow(row_opt[ipoint]); @@ -474,33 +474,33 @@ int main(int argc, char *argv[]) if(z_opt.size()>ipoint){ // assert(x_opt.size()>ipoint); // assert(y_opt.size()>ipoint); - gsl::vector thePosGCP(3); + arma::vec thePosGCP(3); if(x_opt.size()>ipoint&&y_opt.size()>ipoint){ if(gcprad_opt[0]){ - thePosGCP[0]=theModel.rad2deg(x_opt[ipoint]); - thePosGCP[1]=theModel.rad2deg(y_opt[ipoint]); + thePosGCP(0)=theModel.rad2deg(x_opt[ipoint]); + thePosGCP(1)=theModel.rad2deg(y_opt[ipoint]); } else{ - thePosGCP[0]=x_opt[ipoint]; - thePosGCP[1]=y_opt[ipoint]; + thePosGCP(0)=x_opt[ipoint]; + thePosGCP(1)=y_opt[ipoint]; } } - thePosGCP[2]=z_opt[ipoint]; + thePosGCP(2)=z_opt[ipoint]; if(s_srs_opt[0]!=4326) - poCT->Transform(1,&(thePosGCP[0]),&(thePosGCP[1]),&(thePosGCP[2])); + poCT->Transform(1,&(thePosGCP(0)),&(thePosGCP(1)),&(thePosGCP(2))); theDataModel.pushPosGCP(thePosGCP); } assert(yl_opt.size()>ipoint); assert(zl_opt.size()>ipoint); if(pplrad_opt[0]){ - thePosPlatform[0]=theModel.rad2deg(xl_opt[ipoint]); - thePosPlatform[1]=theModel.rad2deg(yl_opt[ipoint]); + thePosPlatform(0)=theModel.rad2deg(xl_opt[ipoint]); + thePosPlatform(1)=theModel.rad2deg(yl_opt[ipoint]); } else{ - thePosPlatform[0]=xl_opt[ipoint]; - thePosPlatform[1]=yl_opt[ipoint]; + thePosPlatform(0)=xl_opt[ipoint]; + thePosPlatform(1)=yl_opt[ipoint]; } - thePosPlatform[2]=zl_opt[ipoint]; + thePosPlatform(2)=zl_opt[ipoint]; // if(s_srs_opt[0]!=4326) // poCT->Transform(1,&(thePosPlatform[0]),&(thePosPlatform[1]),&(thePosPlatform[2])); theDataModel.pushPosPlatform(thePosPlatform); @@ -509,14 +509,14 @@ int main(int argc, char *argv[]) assert(pitch_opt.size()>ipoint); assert(yaw_opt.size()>ipoint); if(aplrad_opt[0]){ - theAttPlatform[0]=theModel.rad2deg(roll_opt[ipoint]); - theAttPlatform[1]=theModel.rad2deg(pitch_opt[ipoint]); - theAttPlatform[2]=theModel.rad2deg(yaw_opt[ipoint]); + theAttPlatform(0)=theModel.rad2deg(roll_opt[ipoint]); + theAttPlatform(1)=theModel.rad2deg(pitch_opt[ipoint]); + theAttPlatform(2)=theModel.rad2deg(yaw_opt[ipoint]); } else{ - theAttPlatform[0]=roll_opt[ipoint]; - theAttPlatform[1]=pitch_opt[ipoint]; - theAttPlatform[2]=yaw_opt[ipoint]; + theAttPlatform(0)=roll_opt[ipoint]; + theAttPlatform(1)=pitch_opt[ipoint]; + theAttPlatform(2)=yaw_opt[ipoint]; } theDataModel.pushAttPlatform(theAttPlatform); } @@ -527,10 +527,10 @@ int main(int argc, char *argv[]) for(int ipoint=0;ipoint<vID.size();++ipoint){ if(verbose_opt[0]>1) std::cout << "point: " << ipoint << std::endl; - gsl::vector pos_platform=theDataModel.getPosPlatform(ipoint); - gsl::vector att_platform=theDataModel.getAttPlatform(ipoint); - gsl::vector pos_calc=theDataModel.getPos(ipoint); - gsl::vector pos_gcp=theDataModel.getPosGCP(ipoint); + arma::vec pos_platform=theDataModel.getPosPlatform(ipoint); + arma::vec att_platform=theDataModel.getAttPlatform(ipoint); + arma::vec pos_calc=theDataModel.getPos(ipoint); + arma::vec pos_gcp=theDataModel.getPosGCP(ipoint); ostringstream gcpss; double e=theDataModel.getModel().getDistGeo(pos_gcp,pos_calc); if(e>=theDataModel.getThreshold()){ @@ -608,11 +608,11 @@ int main(int argc, char *argv[]) for(int ipoint=0;ipoint<vID.size();++ipoint){ if(verbose_opt[0]>1) std::cout << "point: " << ipoint << std::endl; - gsl::vector pos_platform=theDataModel.getPosPlatform(ipoint); - gsl::vector att_platform=theDataModel.getAttPlatform(ipoint); + arma::vec pos_platform=theDataModel.getPosPlatform(ipoint); + arma::vec att_platform=theDataModel.getAttPlatform(ipoint); - gsl::vector pos_calc; - gsl::vector pos_gcp; + arma::vec pos_calc; + arma::vec pos_gcp; if(input_opt.size()||x_opt.size()) pos_gcp=theDataModel.getPosGCP(ipoint); else @@ -628,20 +628,20 @@ int main(int argc, char *argv[]) ostringstream gcpss; gcpss.precision(12); - poCT->Transform(1,&(pos_gcp[0]),&(pos_gcp[1]),&(pos_gcp[2])); - gcpss << pos_gcp[0] << " " << pos_gcp[1] << " " << pos_gcp[2] << " "; + poCT->Transform(1,&(pos_gcp(0)),&(pos_gcp(1)),&(pos_gcp(2))); + gcpss << pos_gcp(0) << " " << pos_gcp(1) << " " << pos_gcp(2) << " "; gcpss << errorv.back() << " "; - poCT->Transform(1,&(pos_calc[0]),&(pos_calc[1]),&(pos_calc[2])); + poCT->Transform(1,&(pos_calc(0)),&(pos_calc(1)),&(pos_calc(2))); if(output_opt.size()){ - outputStream << setprecision(12) << vID[ipoint] << " " << theDataModel.getRow(ipoint) << " " << theDataModel.getCol(ipoint) << " " << pos_calc[0] << " " << pos_calc[1] << " " << pos_calc[2] << " " << gcpss.str(); + outputStream << setprecision(12) << vID[ipoint] << " " << theDataModel.getRow(ipoint) << " " << theDataModel.getCol(ipoint) << " " << pos_calc(0) << " " << pos_calc(1) << " " << pos_calc(2) << " " << gcpss.str(); if(getzenith_opt[0]) outputStream << " " << theDataModel.getModel().getZenith(att_platform,theDataModel.getRow(ipoint),theDataModel.getCol(ipoint)); outputStream << std::endl; } else{ - std::cout << setprecision(12) << vID[ipoint] << " " << theDataModel.getRow(ipoint) << " " << theDataModel.getCol(ipoint) << " " << pos_calc[0] << " " << pos_calc[1] << " " << pos_calc[2] << " " << gcpss.str(); + std::cout << setprecision(12) << vID[ipoint] << " " << theDataModel.getRow(ipoint) << " " << theDataModel.getCol(ipoint) << " " << pos_calc(0) << " " << pos_calc(1) << " " << pos_calc(2) << " " << gcpss.str(); if(getzenith_opt[0]) std::cout << " " << theDataModel.getModel().getZenith(att_platform,theDataModel.getRow(ipoint),theDataModel.getCol(ipoint)); std::cout << std::endl; diff --git a/src/apps/pksensormodel.h b/src/apps/pksensormodel.h index fd13ecb..09a80f4 100644 --- a/src/apps/pksensormodel.h +++ b/src/apps/pksensormodel.h @@ -20,7 +20,8 @@ along with pktools. If not, see <http://www.gnu.org/licenses/>. #ifndef _PKSENSORMODEL_H_ #define _PKSENSORMODEL_H_ #include <vector> -#include <gslwrap/matrix_double.h> +// #include <gslwrap/matrix_double.h> +#include <armadillo> #include "models/SensorModel.h" double objFunction(const std::vector<double> &x, std::vector<double> &grad, void *my_func_data); @@ -42,15 +43,15 @@ class DataModel{ m_row.erase(m_row.begin()+index); m_col.erase(m_col.begin()+index); }; - int pushAttPlatform(const gsl::vector& atp){m_attPlatform.push_back(atp); return m_attPlatform.size();}; - int pushPosPlatform(const gsl::vector& ppl){m_posPlatform.push_back(ppl); return m_posPlatform.size();}; - int pushPosGCP(const gsl::vector& pgcp){m_posGCP.push_back(pgcp); return m_posGCP.size();}; + int pushAttPlatform(const arma::vec& atp){m_attPlatform.push_back(atp); return m_attPlatform.size();}; + int pushPosPlatform(const arma::vec& ppl){m_posPlatform.push_back(ppl); return m_posPlatform.size();}; + int pushPosGCP(const arma::vec& pgcp){m_posGCP.push_back(pgcp); return m_posGCP.size();}; int pushRow(int r){m_row.push_back(r); return m_row.size();}; int pushCol(int c){m_col.push_back(c); return m_col.size();}; - gsl::vector getPosPlatform(int index) const{assert(index>=0);assert(index<m_posPlatform.size());return(m_posPlatform[index]);}; - gsl::vector getAttPlatform(int index) const{assert(index>=0);assert(index<m_attPlatform.size());return(m_attPlatform[index]);}; - gsl::vector getPosGCP(int index) const{assert(index>=0);assert(index<m_posGCP.size());return(m_posGCP[index]);}; - gsl::vector getPos(int index) const{ + arma::vec getPosPlatform(int index) const{assert(index>=0);assert(index<m_posPlatform.size());return(m_posPlatform[index]);}; + arma::vec getAttPlatform(int index) const{assert(index>=0);assert(index<m_attPlatform.size());return(m_attPlatform[index]);}; + arma::vec getPosGCP(int index) const{assert(index>=0);assert(index<m_posGCP.size());return(m_posGCP[index]);}; + arma::vec getPos(int index) const{ assert(index>=0); assert(index<m_posPlatform.size()); assert(index<m_attPlatform.size()); @@ -63,16 +64,16 @@ class DataModel{ int getRow(int index) const{assert(index>=0);assert(index<m_row.size());return(m_row[index]);}; int getCol(int index) const{assert(index>=0);assert(index<m_col.size());return(m_col[index]);}; double getHeight(int index) const{assert(index>=0);assert(index<m_posGCP.size());return(m_posGCP[index][2]);}; - void setBoresightAtt(const gsl::vector& bc_att){ + void setBoresightAtt(const arma::vec& bc_att){ m_model.setBoresightAtt(bc_att); // for(int index=0;index<m_attPlatform.size();++index) // m_attPlatform[index]+=bc_att; }; private: SensorModel::SensorModel m_model; - vector<gsl::vector> m_posPlatform; - vector<gsl::vector> m_posGCP; - vector<gsl::vector> m_attPlatform; + vector<arma::vec> m_posPlatform; + vector<arma::vec> m_posGCP; + vector<arma::vec> m_attPlatform; vector<int> m_row; vector<int> m_col; double m_threshold; diff --git a/src/models/SensorModel.h b/src/models/SensorModel.h index 6bebf9a..77209e2 100644 --- a/src/models/SensorModel.h +++ b/src/models/SensorModel.h @@ -23,16 +23,7 @@ along with pktools. If not, see <http://www.gnu.org/licenses/>. #include <math.h> #include <vector> #include <iostream> -// #include <gsl/gsl_matrix.h> -// #include <gsl/gsl_vector.h> -// #include <gsl/gsl_permutation.h> -// #include <gsl/gsl_linalg.h> -// #include <gsl/gsl_blas.h> -// #include <gsl/gsl_rng.h> -// #include <gsl/gsl_randist.h> -#include <gslwrap/vector_double.h> -#include <gslwrap/matrix_double.h> -#include <gslwrap/matrix_vector_operators.h> +#include <armadillo> #include "ogr_spatialref.h" #ifndef PI @@ -50,8 +41,8 @@ namespace SensorModel class SensorModel{ public: SensorModel(void){ - m_bcpos.resize(3); - m_bcatt.resize(3); + m_bcpos.set_size(3); + m_bcatt.set_size(3); m_bcpos[0]=0; m_bcpos[1]=0; m_bcpos[2]=0; @@ -60,8 +51,8 @@ namespace SensorModel m_bcatt[2]=0; }; SensorModel(int theModel) : m_model(theModel){ - m_bcpos.resize(3); - m_bcatt.resize(3); + m_bcpos.set_size(3); + m_bcatt.set_size(3); m_bcpos[0]=0; m_bcpos[1]=0; m_bcpos[2]=0; @@ -84,117 +75,117 @@ namespace SensorModel double getF() const {return m_fc;}; void setPPx(double ppx){m_ppx=ppx;}; void setPPy(double ppy){m_ppy=ppy;}; - void setBoresightPos(const gsl::vector& bcpos){m_bcpos=bcpos;}; - void setBoresightAtt(const gsl::vector& bcatt){m_bcatt=bcatt;}; - gsl::vector getBoresightPos() const{return m_bcpos;}; - gsl::vector getBoresightAtt() const{return m_bcatt;}; + void setBoresightPos(const arma::vec& bcpos){m_bcpos=bcpos;}; + void setBoresightAtt(const arma::vec& bcatt){m_bcatt=bcatt;}; + arma::vec getBoresightPos() const{return m_bcpos;}; + arma::vec getBoresightAtt() const{return m_bcatt;}; void setPolynome(const std::vector<double>& polynome) {m_polynome=polynome;}; void setDatum(const std::string& theDatum="WGS84"){m_spatialref.SetWellKnownGeogCS(theDatum.c_str());}; - double getZenith(const gsl::vector& att_platform, int row, int column) const{ - gsl::vector normallevel(3); - gsl::vector normalplatform(3); + double getZenith(const arma::vec& att_platform, int row, int column) const{ + arma::vec normallevel(3); + arma::vec normalplatform(3); - gsl::vector apl_deg=att_platform; + arma::vec apl_deg=att_platform; apl_deg+=m_bcatt; - gsl::vector apl_rad(3); - apl_rad[0]=deg2rad(apl_deg[0]);//roll - apl_rad[1]=deg2rad(apl_deg[1]);//pitch - apl_rad[2]=deg2rad(apl_deg[2]);//yaw + arma::vec apl_rad(3); + apl_rad(0)=deg2rad(apl_deg(0));//roll + apl_rad(1)=deg2rad(apl_deg(1));//pitch + apl_rad(2)=deg2rad(apl_deg(2));//yaw if(getModel()==PUSHBROOM){ - gsl::vector scanAngle(2); + arma::vec scanAngle(2); scanAngle=scanAngle_PB(row,column); - apl_rad[0]+=scanAngle[1]; - apl_rad[1]+=scanAngle[0]; + apl_rad(0)+=scanAngle(1); + apl_rad(1)+=scanAngle(0); } else if(getModel()==WHISKBROOM){ - apl_rad[0]+=0; - apl_rad[1]=tan(scanAngle_WB(column)); + apl_rad(0)+=0; + apl_rad(1)=tan(scanAngle_WB(column)); } else//not implemented? assert(0); - normallevel[0]=0; - normallevel[1]=0; - normallevel[2]=-1; - gsl::matrix rotM(3,3); - rotM=getRz(apl_rad[2])*getRy(apl_rad[1])*getRx(apl_rad[0]); + normallevel(0)=0; + normallevel(1)=0; + normallevel(2)=-1; + arma::Mat<double> rotM(3,3); + rotM=getRz(apl_rad(2))*getRy(apl_rad(1))*getRx(apl_rad(0)); normalplatform=rotM*normallevel; - return rad2deg(acos(-normalplatform[2])); + return rad2deg(acos(-normalplatform(2))); }; - gsl::vector getPos(const gsl::vector& pos_platform, const gsl::vector& att_platform, int row, int column, double elevation) const{ - gsl::vector thePosition(3); - gsl::vector pos_ellips(3); - gsl::vector ppl_deg=pos_platform; - gsl::vector apl_deg=att_platform; + arma::vec getPos(const arma::vec& pos_platform, const arma::vec& att_platform, int row, int column, double elevation) const{ + arma::vec thePosition(3); + arma::vec pos_ellips(3); + arma::vec ppl_deg=pos_platform; + arma::vec apl_deg=att_platform; ppl_deg+=m_bcpos; apl_deg+=m_bcatt; - gsl::vector ppl_rad(3); - gsl::vector apl_rad(3); - ppl_rad[0]=deg2rad(ppl_deg[0]); - ppl_rad[1]=deg2rad(ppl_deg[1]); - ppl_rad[2]=ppl_deg[2];//add geoid elevation if necessary... - apl_rad[0]=deg2rad(apl_deg[0]);//roll - apl_rad[1]=deg2rad(apl_deg[1]);//pitch - apl_rad[2]=deg2rad(apl_deg[2]);//yaw - gsl::vector pos_ecef=pECEF(ppl_rad,apl_rad,row,column); + arma::vec ppl_rad(3); + arma::vec apl_rad(3); + ppl_rad(0)=deg2rad(ppl_deg(0)); + ppl_rad(1)=deg2rad(ppl_deg(1)); + ppl_rad(2)=ppl_deg(2);//add geoid elevation if necessary... + apl_rad(0)=deg2rad(apl_deg(0));//roll + apl_rad(1)=deg2rad(apl_deg(1));//pitch + apl_rad(2)=deg2rad(apl_deg(2));//yaw + arma::vec pos_ecef=pECEF(ppl_rad,apl_rad,row,column); pos_ellips=ecef2geo(pos_ecef); - pos_ellips[2]=0; + pos_ellips(2)=0; thePosition=getXYatZ(elevation,ppl_deg,pos_ellips); return thePosition; }; - double getDistGeo(const gsl::vector& pos1, const gsl::vector& pos2) const{ - double lon1=pos1[0]; - double lat1=pos1[1]; - double lon2=pos2[0]; - double lat2=pos2[1]; + double getDistGeo(const arma::vec& pos1, const arma::vec& pos2) const{ + double lon1=pos1(0); + double lat1=pos1(1); + double lon2=pos2(0); + double lat2=pos2(1); double result; //simplified formula (spherical approximation) // result=2*asin(sqrt(pow(sin((lat1-lat2)/2),2) + cos(lat1)*cos(lat2)*pow(sin((lon1-lon2)/2),2))); //using loxodromes - result=(pos1[1]==pos2[1]) ? lox2(deg2rad(pos1[1]),deg2rad(pos1[0]),deg2rad(pos2[0])) : lox1(deg2rad(pos1[1]),deg2rad(pos1[0]),deg2rad(pos2[1]),deg2rad(pos2[0])); + result=(pos1(1)==pos2(1)) ? lox2(deg2rad(pos1(1)),deg2rad(pos1(0)),deg2rad(pos2(0))) : lox1(deg2rad(pos1(1)),deg2rad(pos1(0)),deg2rad(pos2(1)),deg2rad(pos2(0))); return result; }; - gsl::vector geo2ecef(const gsl::vector& pos_geo) const{ - gsl::vector pos_ecef(3); - double nu=getNu(pos_geo[1]); - double f1=(nu+pos_geo[2])*cos(pos_geo[1]); + arma::vec geo2ecef(const arma::vec& pos_geo) const{ + arma::vec pos_ecef(3); + double nu=getNu(pos_geo(1)); + double f1=(nu+pos_geo(2))*cos(pos_geo(1)); double e1=getE1(); - pos_ecef[0]=f1*cos(pos_geo[0]); - pos_ecef[1]=f1*sin(pos_geo[0]); - pos_ecef[2]=(nu*(1-e1*e1)+pos_geo[2])*sin(pos_geo[1]); + pos_ecef(0)=f1*cos(pos_geo(0)); + pos_ecef(1)=f1*sin(pos_geo(0)); + pos_ecef(2)=(nu*(1-e1*e1)+pos_geo(2))*sin(pos_geo(1)); return pos_ecef; }; - gsl::vector ecef2geo(const gsl::vector& pos_ecef) const{ - gsl::vector pos_geo(3); - pos_geo[0]=atan2(pos_ecef[1],pos_ecef[0]); - while(pos_geo[0]>2*PI) - pos_geo[0]-=2*PI; - double r=sqrt(pos_ecef[0]*pos_ecef[0]+pos_ecef[1]*pos_ecef[1]); + arma::vec ecef2geo(const arma::vec& pos_ecef) const{ + arma::vec pos_geo(3); + pos_geo(0)=atan2(pos_ecef(1),pos_ecef(0)); + while(pos_geo(0)>2*PI) + pos_geo(0)-=2*PI; + double r=sqrt(pos_ecef(0)*pos_ecef(0)+pos_ecef(1)*pos_ecef(1)); double f_earth=1.0/m_spatialref.GetInvFlattening(); double e2=f_earth*(2-f_earth); double Ae=m_spatialref.GetSemiMajor(); - pos_geo[1]=atan(pos_ecef[2]/r); + pos_geo(1)=atan(pos_ecef(2)/r); double c=1; double lat=1E+30; int iterations=0; - while(sqrt((pos_geo[1]-lat)*(pos_geo[1]-lat))>1E-15){ + while(sqrt((pos_geo(1)-lat)*(pos_geo(1)-lat))>1E-15){ ++iterations; - lat=pos_geo[1]; + lat=pos_geo(1); c=1.0/sqrt(1-e2*sin(lat)*sin(lat)); - pos_geo[1]=atan((pos_ecef[2]+Ae*c*e2*sin(lat))/r); + pos_geo(1)=atan((pos_ecef(2)+Ae*c*e2*sin(lat))/r); } - pos_geo[2]=r/cos(pos_geo[1])-Ae*c; - pos_geo[0]=rad2deg(pos_geo[0]); - pos_geo[1]=rad2deg(pos_geo[1]); + pos_geo(2)=r/cos(pos_geo(1))-Ae*c; + pos_geo(0)=rad2deg(pos_geo(0)); + pos_geo(1)=rad2deg(pos_geo(1)); return pos_geo; }; private: //line function; In this function, point 0 is the platform position, point 1 is the ellipsoid position - gsl::vector getXYatZ (double inZ, gsl::vector p0, gsl::vector p1) const{ - gsl::vector posatz(3); - posatz[0]=p0[0]+(p0[0]-p1[0])/(p0[2]-p1[2])*(inZ-p0[2]); - posatz[1]=p0[1]+(p0[1]-p1[1])/(p0[2]-p1[2])*(inZ-p0[2]); - posatz[2]=inZ; + arma::vec getXYatZ (double inZ, arma::vec p0, arma::vec p1) const{ + arma::vec posatz(3); + posatz(0)=p0(0)+(p0(0)-p1(0))/(p0(2)-p1(2))*(inZ-p0(2)); + posatz(1)=p0(1)+(p0(1)-p1(1))/(p0(2)-p1(2))*(inZ-p0(2)); + posatz(2)=inZ; return posatz; }; //get First eccentricity of the Earth's surface @@ -228,68 +219,66 @@ namespace SensorModel poCT->Transform(1,&lon,&lat); }; - gsl::vector pECEF(const gsl::vector& pos, const gsl::vector& attitude, int row, int column) const{ - gsl::vector A(3); - gsl::matrix B(3,3); - B.set_element(0,0,-sin(pos[1])*cos(pos[0])); - B.set_element(0,1,-sin(pos[0])); - B.set_element(0,2,-cos(pos[1])*cos(pos[0])); - B.set_element(1,0,-sin(pos[1])*sin(pos[0])); - B.set_element(1,1,cos(pos[0])); - B.set_element(1,2,-cos(pos[1])*sin(pos[0])); - B.set_element(2,0,cos(pos[1])); - B.set_element(2,1,0); - B.set_element(2,2,-sin(pos[1])); + arma::vec pECEF(const arma::vec& pos, const arma::vec& attitude, int row, int column) const{ + arma::vec A(3); + arma::Mat<double> B(3,3); + B(0,0)=-sin(pos(1))*cos(pos(0)); + B(0,1)=-sin(pos(0)); + B(0,2)=-cos(pos(1))*cos(pos(0)); + B(1,0)=-sin(pos(1))*sin(pos(0)); + B(1,1)=cos(pos(0)); + B(1,2)=-cos(pos(1))*sin(pos(0)); + B(2,0)=cos(pos(1)); + B(2,1)=0; + B(2,2)=-sin(pos(1)); A=geo2ecef(pos); - gsl::vector C(3); - A+=B*SV(row,column,attitude,pos[2]); + arma::vec C(3); + A+=B*SV(row,column,attitude,pos(2)); return A; }; - gsl::vector getIV(int row, int column) const{ - gsl::vector iv(3); - iv[0]=0; - gsl::vector scanAngle(2); + arma::vec getIV(int row, int column) const{ + arma::vec iv(3); + iv(0)=0; + arma::vec scanAngle(2); if(getModel()==PUSHBROOM){ scanAngle=scanAngle_PB(row,column); - iv[0]=scanAngle[1]; - iv[1]=tan(scanAngle[0]); + iv(0)=scanAngle(1); + iv(1)=tan(scanAngle(0)); } else if(getModel()==WHISKBROOM){ - iv[0]=0; - iv[1]=tan(scanAngle_WB(column)); + iv(0)=0; + iv(1)=tan(scanAngle_WB(column)); } else//not implemented? assert(0); - iv[2]=-1.0; + iv(2)=-1.0; return iv; }; double scanAngle_WB(int column) const{ return (-m_fov/2.0+column*(m_fov/(m_ncol-1))); }; - gsl::vector SV(int row, int column, const gsl::vector& attitude, double z) const{ - gsl::vector rv(3); - gsl::matrix rotM(3,3); - rotM=getRz(attitude[2])*getRy(attitude[1])*getRx(attitude[0]); + arma::vec SV(int row, int column, const arma::vec& attitude, double z) const{ + arma::vec rv(3); + arma::Mat<double> rotM(3,3); + rotM=getRz(attitude(2))*getRy(attitude(1))*getRx(attitude(0)); rv=rotM*getIV(row,column); - double height=rv[2]; - // gsl_blas_dscal(z/height,&rv); - // return rv; - rv*=z/height; + double height=rv(2); + rv=rv*z/height; return rv; }; - gsl::vector scanAngle_PB(int row, int column) const{ - gsl::vector alpha(2); + arma::vec scanAngle_PB(int row, int column) const{ + arma::vec alpha(2); double r=corr_along(column); double theAx=(column<m_ppx) ? m_dx*0.000001*(m_ppx-(column+corr_across(column))) : m_dx*0.000001*(m_ppx-(column-corr_across(column))); double theAy=m_dy*0.000001*(m_ppy-r); - alpha[0]=(getModel()==FRAME) ? atan(theAx/getF()) : atan(-theAx/getF()); - alpha[1]=atan(-theAy/getF()); + alpha(0)=(getModel()==FRAME) ? atan(theAx/getF()) : atan(-theAx/getF()); + alpha(1)=atan(-theAy/getF()); return alpha; }; //geocentric coordinate system is right-handed, orthogogal Cartesian system with its origin at the centre of the earth. The direct Helmert transformation from lat/lon and ellipsoid height to geocentric x,y,z is returned in posX, posY, posZ //Euclidian distance between points on sphere //Euclidian distance between two vectors - double getDist(const gsl::vector& v1, const gsl::vector& v2) const{gsl::vector dv=v1;dv-=v2;return dv.norm2();}; + double getDist(const arma::vec& v1, const arma::vec& v2) const{arma::vec dv=v1;dv-=v2;return arma::norm(dv,2);}; double lox1(double lat1, double lon1, double lat2, double lon2) const{ double result=(M(lat2)-M(lat1))/cos(Az(lat1,lon1,lat2,lon2)); return result; @@ -348,43 +337,44 @@ namespace SensorModel } return result; }; - gsl::matrix getRx(double theta) const{ - gsl::matrix Rx(3,3); - Rx.set_element(0,0,1); - Rx.set_element(0,1,0); - Rx.set_element(0,2,0); - Rx.set_element(1,0,0); - Rx.set_element(1,1,cos(theta)); - Rx.set_element(1,2,-sin(theta)); - Rx.set_element(2,0,0); - Rx.set_element(2,1,sin(theta)); - Rx.set_element(2,2,cos(theta)); + arma::Mat<double> getRx(double theta) const{ + arma::Mat<double> Rx(3,3); + Rx(0,0)=1; + Rx(0,1)=0; + Rx(0,2)=0; + Rx(1,0)=0; + Rx(1,1)=cos(theta); + Rx(1,2)=-sin(theta); + Rx(2,0)=0; + Rx(2,1)=sin(theta); + Rx(2,2)=cos(theta); return Rx; }; - gsl::matrix getRy(double theta) const{ - gsl::matrix Ry(3,3); - Ry.set_element(0,0,cos(theta)); - Ry.set_element(0,1,0); - Ry.set_element(0,2,sin(theta)); - Ry.set_element(1,0,0); - Ry.set_element(1,1,1); - Ry.set_element(1,2,0); - Ry.set_element(2,0,-sin(theta)); - Ry.set_element(2,1,0); - Ry.set_element(2,2,cos(theta)); + arma::Mat<double> getRy(double theta) const{ + arma::Mat<double> Ry(3,3); + Ry(0,0)=cos(theta); + Ry(0,1)=0; + Ry(0,2)=sin(theta); + Ry(1,0)=0; + Ry(1,1)=1; + Ry(1,2)=0; + Ry(2,0)=-sin(theta); + Ry(2,1)=0; + Ry(2,2)=cos(theta); return Ry; }; - gsl::matrix getRz(double theta) const{ - gsl::matrix Rz(3,3); - Rz.set_element(0,0,cos(theta)); - Rz.set_element(0,1,-sin(theta)); - Rz.set_element(0,2,0); - Rz.set_element(1,0,sin(theta)); - Rz.set_element(1,1,cos(theta)); - Rz.set_element(1,2,0); - Rz.set_element(2,0,0); - Rz.set_element(2,1,0); - Rz.set_element(2,2,1); + + arma::Mat<double> getRz(double theta) const{ + arma::Mat<double> Rz(3,3); + Rz(0,0)=cos(theta); + Rz(0,1)=-sin(theta); + Rz(0,2)=0; + Rz(1,0)=sin(theta); + Rz(1,1)=cos(theta); + Rz(1,2)=0; + Rz(2,0)=0; + Rz(2,1)=0; + Rz(2,2)=1; return Rz; }; int m_model; @@ -401,8 +391,8 @@ namespace SensorModel int m_nrow; double m_dx; double m_dy; - gsl::vector m_bcpos; - gsl::vector m_bcatt; + arma::vec m_bcpos; + arma::vec m_bcatt; std::vector<double> m_polynome; }; } -- Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-grass/pktools.git _______________________________________________ Pkg-grass-devel mailing list Pkg-grass-devel@lists.alioth.debian.org http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/pkg-grass-devel