[pktools] 253/375: implementing vito dsm2dtm
Bas Couwenberg
sebastic at xs4all.nl
Wed Dec 3 21:54:19 UTC 2014
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 5494440653ecbce631b45b1067c595f40857c7ce
Author: Pieter Kempeneers <kempenep at gmail.com>
Date: Wed Apr 30 23:00:28 2014 +0200
implementing vito dsm2dtm
---
src/algorithms/Filter2d.h | 24 ++++++----
src/apps/pkfilterdem.cc | 111 ++++++++++++++++++++++++++++++++++++++++++----
src/apps/pklas2img.cc | 12 ++---
src/base/Vector2d.h | 21 +++++++++
4 files changed, 146 insertions(+), 22 deletions(-)
diff --git a/src/algorithms/Filter2d.h b/src/algorithms/Filter2d.h
index bc33cfc..54aeee0 100644
--- a/src/algorithms/Filter2d.h
+++ b/src/algorithms/Filter2d.h
@@ -825,12 +825,14 @@ template<class T> unsigned long int Filter2d::morphology(const Vector2d<T>& inpu
++nmasked;
}
}
- if(nmasked>=nlimit){
+ if(nmasked<nlimit){
++nchange;
//reset pixel in outputMask
outputMask[y][x]=0;
+ }
+ else{
//reset pixel height in tmpDSM
- inBuffer[indexJ][indexI]=stat.mymin(neighbors);
+ inBuffer[(dimY-1)/2][x]=stat.mymin(neighbors);
}
}
progress=(1.0+y);
@@ -914,12 +916,14 @@ template<class T> unsigned long int Filter2d::morphology(const Vector2d<T>& inpu
++nmasked;
}
}
- if(nmasked>=nlimit){
+ if(nmasked<nlimit){
++nchange;
//reset pixel in outputMask
outputMask[y][x]=0;
+ }
+ else{
//reset pixel height in tmpDSM
- inBuffer[indexJ][indexI]=stat.mymin(neighbors);
+ inBuffer[(dimY-1)/2][x]=stat.mymin(neighbors);
}
}
progress=(1.0+y);
@@ -1003,12 +1007,14 @@ template<class T> unsigned long int Filter2d::morphology(const Vector2d<T>& inpu
++nmasked;
}
}
- if(nmasked>=nlimit){
+ if(nmasked<nlimit){
++nchange;
//reset pixel in outputMask
outputMask[y][x]=0;
+ }
+ else{
//reset pixel height in tmpDSM
- inBuffer[indexJ][indexI]=stat.mymin(neighbors);
+ inBuffer[(dimY-1)/2][x]=stat.mymin(neighbors);
}
}
progress=(1.0+y);
@@ -1092,12 +1098,14 @@ template<class T> unsigned long int Filter2d::morphology(const Vector2d<T>& inpu
++nmasked;
}
}
- if(nmasked>=nlimit){
+ if(nmasked<nlimit){
++nchange;
//reset pixel in outputMask
outputMask[y][x]=0;
+ }
+ else{
//reset pixel height in tmpDSM
- inBuffer[indexJ][indexI]=stat.mymin(neighbors);
+ inBuffer[(dimY-1)/2][x]=stat.mymin(neighbors);
}
}
progress=(1.0+y);
diff --git a/src/apps/pkfilterdem.cc b/src/apps/pkfilterdem.cc
index 3238264..b96e06c 100644
--- a/src/apps/pkfilterdem.cc
+++ b/src/apps/pkfilterdem.cc
@@ -161,32 +161,127 @@ int main(int argc,char **argv) {
unsigned long int nchange=1;
if(postFilter_opt[0]=="vito"){
//todo: fill empty pixels
- hThreshold_opt.resize(3);
- hThreshold_opt[0]=0.7;
- hThreshold_opt[1]=0.3;
- hThreshold_opt[2]=0.1;
+ // hThreshold_opt.resize(3);
+ // hThreshold_opt[0]=0.7;
+ // hThreshold_opt[1]=0.3;
+ // hThreshold_opt[2]=0.1;
vector<int> nlimit(3);
nlimit[0]=2;
nlimit[1]=3;
nlimit[2]=4;
+ nlimit[2]=2;
//init finalMask
for(int irow=0;irow<tmpData.nRows();++irow)
for(int icol=0;icol<tmpData.nCols();++icol)
tmpData[irow][icol]=1;
- for(int iheight=0;iheight=hThreshold_opt.size();++iheight){
- //todo: init tmpMask to 1
+ for(int iheight=0;iheight<hThreshold_opt.size();++iheight){
+ if(verbose_opt[0])
+ cout << "height: " << hThreshold_opt[iheight] << endl;
//todo:replace with binary mask (or short) -> adapt template with T1,T2 in Filter2d
Vector2d<double> tmpMask(input.nrOfRow(),input.nrOfCol());
for(int irow=0;irow<tmpMask.nRows();++irow)
for(int icol=0;icol<tmpMask.nCols();++icol)
tmpMask[irow][icol]=1;//1=surface, 0=terrain
+ if(verbose_opt[0])
+ cout << "filtering NWSE" << endl;
theFilter.dsm2dtm_nwse(inputData,tmpMask,hThreshold_opt[iheight],nlimit[iheight],dim_opt[0]);
+
+ // //from here
+
+ // Vector2d<double> tmpDSM(inputData);
+ // double noDataValue=0;
+
+ // unsigned long int nchange=0;
+ // int dimX=dim_opt[0];
+ // int dimY=dim_opt[0];
+ // assert(dimX);
+ // assert(dimY);
+ // statfactory::StatFactory stat;
+ // Vector2d<double> inBuffer(dimY,inputData.nCols());
+ // // if(tmpMask.size()!=inputData.nRows())
+ // // tmpMask.resize(inputData.nRows());
+ // int indexI=0;
+ // int indexJ=0;
+ // //initialize last half of inBuffer
+ // for(int j=-(dimY-1)/2;j<=dimY/2;++j){
+ // for(int i=0;i<inputData.nCols();++i)
+ // inBuffer[indexJ][i]=tmpDSM[abs(j)][i];
+ // ++indexJ;
+ // }
+ // for(int y=0;y<tmpDSM.nRows();++y){
+ // if(y){//inBuffer already initialized for y=0
+ // //erase first line from inBuffer
+ // inBuffer.erase(inBuffer.begin());
+ // //read extra line and push back to inBuffer if not out of bounds
+ // if(y+dimY/2<tmpDSM.nRows()){
+ // //allocate buffer
+ // inBuffer.push_back(inBuffer.back());
+ // for(int i=0;i<tmpDSM.nCols();++i)
+ // inBuffer[inBuffer.size()-1][i]=tmpDSM[y+dimY/2][i];
+ // }
+ // else{
+ // int over=y+dimY/2-tmpDSM.nRows();
+ // int index=(inBuffer.size()-1)-over;
+ // assert(index>=0);
+ // assert(index<inBuffer.size());
+ // inBuffer.push_back(inBuffer[index]);
+ // }
+ // }
+ // for(int x=0;x<tmpDSM.nCols();++x){
+ // double centerValue=inBuffer[(dimY-1)/2][x];
+ // short nmasked=0;
+ // std::vector<double> neighbors;
+ // for(int j=-(dimY-1)/2;j<=dimY/2;++j){
+ // for(int i=-(dimX-1)/2;i<=dimX/2;++i){
+ // indexI=x+i;
+ // //check if out of bounds
+ // if(indexI<0)
+ // indexI=-indexI;
+ // else if(indexI>=tmpDSM.nCols())
+ // indexI=tmpDSM.nCols()-i;
+ // if(y+j<0)
+ // indexJ=-j;
+ // else if(y+j>=tmpDSM.nRows())
+ // indexJ=(dimY>2) ? (dimY-1)/2-j : 0;
+ // else
+ // indexJ=(dimY-1)/2+j;
+ // double difference=(centerValue-inBuffer[indexJ][indexI]);
+ // if(i||j)//skip centerValue
+ // neighbors.push_back(inBuffer[indexJ][indexI]);
+ // if(difference>hThreshold_opt[iheight])
+ // ++nmasked;
+ // }
+ // }
+ // if(nmasked<nlimit[iheight]){
+ // ++nchange;
+ // //reset pixel in tmpMask
+ // tmpMask[y][x]=0;
+ // }
+ // else{
+ // //reset pixel height in tmpDSM
+ // inBuffer[(dimY-1)/2][x]=stat.mymin(neighbors);
+ // }
+ // }
+ // }
+ //to here
+
+ tmpData.setMask(tmpMask,0,0);
+ if(verbose_opt[0])
+ cout << "filtering NESW" << endl;
theFilter.dsm2dtm_nesw(inputData,tmpMask,hThreshold_opt[iheight],nlimit[iheight],dim_opt[0]);
+ tmpData.setMask(tmpMask,0,0);
+ if(verbose_opt[0])
+ cout << "filtering SENW" << endl;
theFilter.dsm2dtm_senw(inputData,tmpMask,hThreshold_opt[iheight],nlimit[iheight],dim_opt[0]);
+ tmpData.setMask(tmpMask,0,0);
+ if(verbose_opt[0])
+ cout << "filtering SWNE" << endl;
theFilter.dsm2dtm_swne(inputData,tmpMask,hThreshold_opt[iheight],nlimit[iheight],dim_opt[0]);
- //todo: set tmpMask to finalMask
+ // set tmpMask to finalMask
+ tmpData.setMask(tmpMask,0,0);
}
- //todo: fill empty pixels
+ outputData=tmpData;
+ //outputData.setMask(tmpData,1,0);
}
else if(postFilter_opt[0]=="etew_min"){
//Elevation Threshold with Expand Window (ETEW) Filter (p.73 from Airborne LIDAR Data Processing and Analysis Tools ALDPAT 1.0)
diff --git a/src/apps/pklas2img.cc b/src/apps/pklas2img.cc
index d0b3404..dad94b0 100644
--- a/src/apps/pklas2img.cc
+++ b/src/apps/pklas2img.cc
@@ -40,7 +40,7 @@ int main(int argc,char **argv) {
Optionpk<unsigned short> returns_opt("ret", "ret", "number(s) of returns to include");
Optionpk<unsigned short> classes_opt("class", "class", "classes to keep: 0 (created, never classified), 1 (unclassified), 2 (ground), 3 (low vegetation), 4 (medium vegetation), 5 (high vegetation), 6 (building), 7 (low point, noise), 8 (model key-point), 9 (water), 10 (reserved), 11 (reserved), 12 (overlap)");
Optionpk<string> composite_opt("comp", "comp", "composite for multiple points in cell (min, max, median, mean, sum, first, last, profile (percentile height values), number (point density)). Last: overwrite cells with latest point", "last");
- Optionpk<string> filter_opt("fir", "filter", "filter las points (last,single,multiple,all).", "all");
+ Optionpk<string> filter_opt("fir", "filter", "filter las points (first,last,single,multiple,all).", "all");
// Optionpk<string> postFilter_opt("pf", "pfilter", "post processing filter (etew_min,promorph (progressive morphological filter),bunting (adapted promorph),open,close,none).", "none");
// Optionpk<short> dimx_opt("dimx", "dimx", "Dimension X of postFilter", 3);
// Optionpk<short> dimy_opt("dimy", "dimy", "Dimension Y of postFilter", 3);
@@ -99,7 +99,7 @@ int main(int argc,char **argv) {
std::cout << "short option -h shows basic options only, use long option --help to show all options" << std::endl;
exit(0);//help was invoked, stop processing
}
-
+ //todo: is this needed?
GDALAllRegister();
double dfComplete=0.0;
@@ -294,10 +294,10 @@ int main(int argc,char **argv) {
continue;
if((filter_opt[0]=="multiple")&&(thePoint.GetNumberOfReturns()<2))
continue;
- if(filter_opt[0]=="last"){
- if(thePoint.GetReturnNumber()!=thePoint.GetNumberOfReturns())
- continue;
- }
+ if((filter_opt[0]=="last")&&(thePoint.GetReturnNumber()!=thePoint.GetNumberOfReturns()))
+ continue;
+ if((filter_opt[0]=="first")&&(thePoint.GetReturnNumber()!=1))
+ continue;
double dcol,drow;
outputWriter.geo2image(thePoint.GetX(),thePoint.GetY(),dcol,drow);
int icol=static_cast<int>(dcol);
diff --git a/src/base/Vector2d.h b/src/base/Vector2d.h
index 12ac813..bbbcafa 100644
--- a/src/base/Vector2d.h
+++ b/src/base/Vector2d.h
@@ -49,6 +49,7 @@ public:
void selectCol(int col, T* output) const;
std::vector<T> selectCol(int col);
void selectCols(const std::list<int> &cols, Vector2d<T> &output) const;
+ void setMask(const Vector2d<T> &mask, T msknodata, T nodata=0);
void transpose(Vector2d<T> &output) const{
output.resize(nCols(),nRows());
for(int irow=0;irow<nRows();++irow){
@@ -62,6 +63,7 @@ public:
void scale(const std::vector<double> &scaleVector, const std::vector<double> &offsetVector, Vector2d<T>& scaledOutput);
void scale(const T lbound, const T ubound, std::vector<double> &scaleVector, std::vector<double> &offsetVector, Vector2d<T>& scaledOutput);
Vector2d<T> operator=(const Vector2d<T>& v1);
+ Vector2d<T> operator+=(const Vector2d<T>& v1);
// std::ostream& operator<<(std::ostream& os, const Vector2d<T>& v);
// template<class T> std::ostream& operator<<(std::ostream& os, const Vector2d<T>& v);
template<class T1> friend std::ostream& operator<<(std::ostream & os, const Vector2d<T1>& v);
@@ -99,6 +101,15 @@ template<class T> Vector2d<T> Vector2d<T>::operator=(const Vector2d<T>& v1){
}
}
+template<class T> Vector2d<T> Vector2d<T>::operator+=(const Vector2d<T>& v1){
+ assert(v1.nRows()==nRows());
+ assert(v1.nCols()==nCols());
+ for(int irow=0;irow<nRows();++irow)
+ for(int icol=0;icol<nCols();++icol)
+ (*this)[irow][icol]+=v1[irow][icol];
+ return *this;
+}
+
template<class T> Vector2d<T>::Vector2d(int nrow)
: std::vector< std::vector<T> >(nrow)
{
@@ -192,6 +203,16 @@ template<class T> void Vector2d<T>::selectCols(const std::list<int> &cols)
(*this)[irow].erase(((*this)[irow]).begin()+icol);
}
+template<class T> void Vector2d<T>::setMask(const Vector2d<T> &mask, T msknodata, T nodata)
+{
+ assert(mask.nRows()==nRows());
+ assert(mask.nCols()==nCols());
+ for(int irow=0;irow<this->size();++irow)
+ for(int icol=0;icol<((*this)[irow]).size()-1;++icol)
+ if(mask[irow][icol]==msknodata)
+ (*this)[irow][icol]=nodata;
+}
+
template<class T1> std::ostream& operator<<(std::ostream& os, const Vector2d<T1>& v)
{
for(int irow=0;irow<v.size();++irow){
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/pkg-grass/pktools.git
More information about the Pkg-grass-devel
mailing list