[pktools] 39/375: see Changelog
Bas Couwenberg
sebastic at xs4all.nl
Wed Dec 3 21:53:56 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 e93a8388d46cc553b26f827aded82d1f9bd89e52
Author: Pieter Kempeneers <kempenep at gmail.com>
Date: Fri Jan 25 17:39:34 2013 +0100
see Changelog
---
ChangeLog | 11 +-
src/algorithms/Filter2d.cc | 76 ++++++++---
src/algorithms/Filter2d.h | 281 +++++++++++++++++++++++++++++++++++++--
src/apps/pkclassify_svm.cc | 1 +
src/apps/pkcrop.cc | 4 +-
src/apps/pkfilter.cc | 24 ++--
src/apps/pklas2img.cc | 26 ++--
src/fileclasses/FileReaderLas.cc | 39 ++++++
src/fileclasses/FileReaderLas.h | 12 ++
9 files changed, 422 insertions(+), 52 deletions(-)
diff --git a/ChangeLog b/ChangeLog
index 93b1c13..f907a49 100755
--- a/ChangeLog
+++ b/ChangeLog
@@ -42,14 +42,20 @@ version 2.4
--enable-nlopt (when NLOPT is installed, needed for pksensormodel and pkgetchandelier)
- OptFactory.h (added)
factory class for nlopt::opt (selecting algorithm via string)
+ - Filter2d
+ support filtering of matrix (doit with Vector2d arguments) (to create DTM according to P. Bunting)
+ - FileReaderLas
+ support class filters
- ImgReaderGdal.cc
in addition to internal setNoData member variable, also support GDALSetNoData
- ImgWriterGdal.cc
in addition to internal setNoData member variable, also support GDALSetNoData
- ImgWriterOGR.h
renamed ascii2shape to ascii2ogr, support csv file in ascii2ogr
- - ImgWriterOGR.cc
+ - ImgWriterOgr.cc
renamed ascii2shape to ascii2ogr, support csv file in ascii2ogr
+ - pkcrop
+ changed default option for color table to empty vector
- pkinfo
support of nodata value when calculating histogram and image statistics
options min and max are now set with -min (--min) and -max (--max)
@@ -58,6 +64,7 @@ version 2.4
- pkfilter
bug solved in update of progress bar
support of standard deviation
+ default empty classes
- pkgetmask
options min and max are now set with -min (--min) and -max (--max)
- pkstatogr
@@ -81,3 +88,5 @@ version 2.4
tool to model pushbroom sensor (with optimization of boresight angles using NLOPT)
- pkascii2ogr
support csv input file
+ - pklas2img
+ support DTM according to Bunting (slightly adapted version of Chang2003, including a median filter after morphological operator)
diff --git a/src/algorithms/Filter2d.cc b/src/algorithms/Filter2d.cc
index 12dd080..e46af52 100644
--- a/src/algorithms/Filter2d.cc
+++ b/src/algorithms/Filter2d.cc
@@ -503,53 +503,85 @@ void Filter2d::Filter2d::doit(const ImgReaderGdal& input, ImgWriterGdal& output,
}
switch(method){
case(MEDIAN):
- outBuffer[x/down]=hist.median(windowBuffer);
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.median(windowBuffer);
break;
case(VAR):{
- outBuffer[x/down]=hist.var(windowBuffer);
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.var(windowBuffer);
break;
}
case(STDEV):{
- outBuffer[x/down]=sqrt(hist.var(windowBuffer));
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=sqrt(hist.var(windowBuffer));
break;
}
case(MEAN):{
- outBuffer[x/down]=hist.mean(windowBuffer);
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.mean(windowBuffer);
break;
}
case(MIN):{
- outBuffer[x/down]=hist.min(windowBuffer);
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.min(windowBuffer);
break;
}
case(ISMIN):{
- outBuffer[x/down]=(hist.min(windowBuffer)==windowBuffer[dimX*dimY/2])? 1:0;
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=(hist.min(windowBuffer)==windowBuffer[dimX*dimY/2])? 1:0;
break;
}
case(MINMAX):{
double min=0;
double max=0;
- hist.minmax(windowBuffer,windowBuffer.begin(),windowBuffer.end(),min,max);
- if(min!=max)
- outBuffer[x/down]=0;
- else
- outBuffer[x/down]=windowBuffer[dimX*dimY/2];//centre pixels
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else{
+ hist.minmax(windowBuffer,windowBuffer.begin(),windowBuffer.end(),min,max);
+ if(min!=max)
+ outBuffer[x/down]=0;
+ else
+ outBuffer[x/down]=windowBuffer[dimX*dimY/2];//centre pixels
+ }
break;
}
case(MAX):{
- outBuffer[x/down]=hist.max(windowBuffer);
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.max(windowBuffer);
break;
}
case(ISMAX):{
- outBuffer[x/down]=(hist.max(windowBuffer)==windowBuffer[dimX*dimY/2])? 1:0;
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=(hist.max(windowBuffer)==windowBuffer[dimX*dimY/2])? 1:0;
break;
}
case(ORDER):{
- double lbound=0;
- double ubound=dimX*dimY;
- double theMin=hist.min(windowBuffer);
- double theMax=hist.max(windowBuffer);
- double scale=(ubound-lbound)/(theMax-theMin);
- outBuffer[x/down]=static_cast<short>(scale*(windowBuffer[dimX*dimY/2]-theMin)+lbound);
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else{
+ double lbound=0;
+ double ubound=dimX*dimY;
+ double theMin=hist.min(windowBuffer);
+ double theMax=hist.max(windowBuffer);
+ double scale=(ubound-lbound)/(theMax-theMin);
+ outBuffer[x/down]=static_cast<short>(scale*(windowBuffer[dimX*dimY/2]-theMin)+lbound);
+ }
break;
}
case(SUM):{
@@ -582,7 +614,7 @@ void Filter2d::Filter2d::doit(const ImgReaderGdal& input, ImgWriterGdal& output,
outBuffer[x/down]+=100.0*occurrence[*(vit++)]/windowBuffer.size();
}
else
- outBuffer[x/down]=0;
+ outBuffer[x/down]=m_noValue;
break;
}
case(MAJORITY):{
@@ -598,7 +630,7 @@ void Filter2d::Filter2d::doit(const ImgReaderGdal& input, ImgWriterGdal& output,
outBuffer[x/down]=inBuffer[dimY/2][x];
}
else
- outBuffer[x/down]=0;
+ outBuffer[x/down]=m_noValue;
break;
}
case(THRESHOLD):{
@@ -611,7 +643,7 @@ void Filter2d::Filter2d::doit(const ImgReaderGdal& input, ImgWriterGdal& output,
}
}
else
- outBuffer[x/down]=0;
+ outBuffer[x/down]=m_noValue;
break;
}
case(MIXED):{
diff --git a/src/algorithms/Filter2d.h b/src/algorithms/Filter2d.h
index 3a76df0..954f1fc 100644
--- a/src/algorithms/Filter2d.h
+++ b/src/algorithms/Filter2d.h
@@ -69,6 +69,7 @@ public:
/* void homogeneousSpatial(const string& inputFilename, const string& outputFilename, int dim, bool disc=false, int noValue=0); */
void doit(const ImgReaderGdal& input, ImgWriterGdal& output, int method, int dim, short down=2, bool disc=false);
void doit(const ImgReaderGdal& input, ImgWriterGdal& output, int method, int dimX, int dimY, short down=1, bool disc=false);
+ template<class T1, class T2> void doit(const Vector2d<T1>& inputVector, Vector2d<T2>& outputVector, int method, int dimX, int dimY, short down=1, bool disc=false);
void median(const string& inputFilename, const string& outputFilename, int dim, bool disc=false);
void var(const string& inputFilename, const string& outputFilename, int dim, bool disc=false);
void morphology(const ImgReaderGdal& input, ImgWriterGdal& output, int method, int dimX, int dimY, bool disc=false, double angle=-190);
@@ -101,13 +102,59 @@ private:
filter(inputVector,outputVector);
}
- template<class T1, class T2> void Filter2d::filter(const Vector2d<T1>& inputVector, Vector2d<T2>& outputVector)
+ template<class T1, class T2> void Filter2d::filter(const Vector2d<T1>& inputVector, Vector2d<T2>& outputVector)
{
- outputVector.resize(inputVector.size());
- int dimX=m_taps[0].size();//horizontal!!!
- int dimY=m_taps.size();//vertical!!!
+ outputVector.resize(inputVector.size());
+ int dimX=m_taps[0].size();//horizontal!!!
+ int dimY=m_taps.size();//vertical!!!
+ Vector2d<T1> inBuffer(dimY);
+ vector<T2> outBuffer(inputVector[0].size());
+ //initialize last half of inBuffer
+ int indexI=0;
+ int indexJ=0;
+ for(int y=0;y<dimY;++y){
+ if(y<dimY/2)
+ continue;//skip first half
+ inBuffer[y]=inputVector[indexJ++];
+ }
+ for(int y=0;y<inputVector.size();++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<inputVector.size())
+ inBuffer.push_back(inputVector[y+dimY/2]);
+ }
+ for(int x=0;x<inputVector[0].size();++x){
+ outBuffer[x]=0;
+ for(int j=-dimY/2;j<(dimY+1)/2;++j){
+ for(int i=-dimX/2;i<(dimX+1)/2;++i){
+ indexI=x+i;
+ indexJ=dimY/2+j;
+ //check if out of bounds
+ if(x<dimX/2)
+ indexI=x+abs(i);
+ else if(x>=inputVector[0].size()-dimX/2)
+ indexI=x-abs(i);
+ if(y<dimY/2)
+ indexJ=dimY/2+abs(j);
+ else if(y>=inputVector.size()-dimY/2)
+ indexJ=dimY/2-abs(j);
+ outBuffer[x]+=(m_taps[dimY/2+j][dimX/2+i]*inBuffer[indexJ][indexI]);
+ }
+ }
+ }
+ //copy outBuffer to outputVector
+ outputVector[y]=outBuffer;
+ }
+ }
+
+template<class T1, class T2> void Filter2d::doit(const Vector2d<T1>& inputVector, Vector2d<T2>& outputVector, int method, int dimX, int dimY, short down, bool disc)
+{
+ Histogram hist;
+ outputVector.resize((inputVector.size()+down-1)/down);
Vector2d<T1> inBuffer(dimY);
- vector<T2> outBuffer(inputVector[0].size());
+ vector<T2> outBuffer((inputVector[0].size()+down-1)/down);
//initialize last half of inBuffer
int indexI=0;
int indexJ=0;
@@ -116,6 +163,11 @@ private:
continue;//skip first half
inBuffer[y]=inputVector[indexJ++];
}
+ const char* pszMessage;
+ void* pProgressArg=NULL;
+ GDALProgressFunc pfnProgress=GDALTermProgress;
+ double progress=0;
+ pfnProgress(progress,pszMessage,pProgressArg);
for(int y=0;y<inputVector.size();++y){
if(y){//inBuffer already initialized for y=0
//erase first line from inBuffer
@@ -123,9 +175,22 @@ private:
//read extra line and push back to inBuffer if not out of bounds
if(y+dimY/2<inputVector.size())
inBuffer.push_back(inputVector[y+dimY/2]);
+ else{
+ int over=y+dimY/2-inputVector.size();
+ int index=(inBuffer.size()-1)-over;
+ assert(index>=0);
+ assert(index<inBuffer.size());
+ inBuffer.push_back(inBuffer[index]);
+ }
}
+ if((y+1+down/2)%down)
+ continue;
for(int x=0;x<inputVector[0].size();++x){
- outBuffer[x]=0;
+ if((x+1+down/2)%down)
+ continue;
+ outBuffer[x/down]=0;
+ vector<double> windowBuffer;
+ map<int,int> occurrence;
for(int j=-dimY/2;j<(dimY+1)/2;++j){
for(int i=-dimX/2;i<(dimX+1)/2;++i){
indexI=x+i;
@@ -139,12 +204,210 @@ private:
indexJ=dimY/2+abs(j);
else if(y>=inputVector.size()-dimY/2)
indexJ=dimY/2-abs(j);
- outBuffer[x]+=(m_taps[dimY/2+j][dimX/2+i]*inBuffer[indexJ][indexI]);
- }
+ bool masked=false;
+ for(int imask=0;imask<m_mask.size();++imask){
+ if(inBuffer[indexJ][indexI]==m_mask[imask]){
+ masked=true;
+ break;
+ }
+ }
+ if(!masked){
+ vector<short>::const_iterator vit=m_class.begin();
+ //todo: test if this works (only add occurrence if within defined classes)!
+ if(!m_class.size())
+ ++occurrence[inBuffer[indexJ][indexI]];
+ else{
+ while(vit!=m_class.end()){
+ if(inBuffer[indexJ][indexI]==*(vit++))
+ ++occurrence[inBuffer[indexJ][indexI]];
+ }
+ }
+ windowBuffer.push_back(inBuffer[indexJ][indexI]);
+ }
+ }
+ }
+ switch(method){
+ case(MEDIAN):
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.median(windowBuffer);
+ break;
+ case(VAR):{
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.var(windowBuffer);
+ break;
+ }
+ case(STDEV):{
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=sqrt(hist.var(windowBuffer));
+ break;
+ }
+ case(MEAN):{
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.mean(windowBuffer);
+ break;
+ }
+ case(MIN):{
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.min(windowBuffer);
+ break;
+ }
+ case(ISMIN):{
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=(hist.min(windowBuffer)==windowBuffer[dimX*dimY/2])? 1:0;
+ break;
+ }
+ case(MINMAX):{
+ double min=0;
+ double max=0;
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else{
+ hist.minmax(windowBuffer,windowBuffer.begin(),windowBuffer.end(),min,max);
+ if(min!=max)
+ outBuffer[x/down]=0;
+ else
+ outBuffer[x/down]=windowBuffer[dimX*dimY/2];//centre pixels
+ }
+ break;
+ }
+ case(MAX):{
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=hist.max(windowBuffer);
+ break;
+ }
+ case(ISMAX):{
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else
+ outBuffer[x/down]=(hist.max(windowBuffer)==windowBuffer[dimX*dimY/2])? 1:0;
+ break;
+ }
+ case(ORDER):{
+ if(windowBuffer.empty())
+ outBuffer[x/down]=m_noValue;
+ else{
+ double lbound=0;
+ double ubound=dimX*dimY;
+ double theMin=hist.min(windowBuffer);
+ double theMax=hist.max(windowBuffer);
+ double scale=(ubound-lbound)/(theMax-theMin);
+ outBuffer[x/down]=static_cast<short>(scale*(windowBuffer[dimX*dimY/2]-theMin)+lbound);
+ }
+ break;
+ }
+ case(SUM):{
+ outBuffer[x/down]=hist.sum(windowBuffer);
+ break;
+ }
+ case(HOMOG):
+ if(occurrence.size()==1)//all values in window must be the same
+ outBuffer[x/down]=inBuffer[dimY/2][x];
+ else//favorize original value in case of ties
+ outBuffer[x/down]=m_noValue;
+ break;
+ case(HETEROG):{
+ for(vector<double>::const_iterator wit=windowBuffer.begin();wit!=windowBuffer.end();++wit){
+ if(wit==windowBuffer.begin()+windowBuffer.size()/2)
+ continue;
+ else if(*wit!=inBuffer[dimY/2][x])
+ outBuffer[x/down]=1;
+ else if(*wit==inBuffer[dimY/2][x]){//todo:wit mag niet central pixel zijn
+ outBuffer[x/down]=m_noValue;
+ break;
+ }
+ }
+ break;
+ }
+ case(DENSITY):{
+ if(windowBuffer.size()){
+ vector<short>::const_iterator vit=m_class.begin();
+ while(vit!=m_class.end())
+ outBuffer[x/down]+=100.0*occurrence[*(vit++)]/windowBuffer.size();
+ }
+ else
+ outBuffer[x/down]=m_noValue;
+ break;
+ }
+ case(MAJORITY):{
+ if(occurrence.size()){
+ map<int,int>::const_iterator maxit=occurrence.begin();
+ for(map<int,int>::const_iterator mit=occurrence.begin();mit!=occurrence.end();++mit){
+ if(mit->second>maxit->second)
+ maxit=mit;
+ }
+ if(occurrence[inBuffer[dimY/2][x]]<maxit->second)//
+ outBuffer[x/down]=maxit->first;
+ else//favorize original value in case of ties
+ outBuffer[x/down]=inBuffer[dimY/2][x];
+ }
+ else
+ outBuffer[x/down]=m_noValue;
+ break;
+ }
+ case(THRESHOLD):{
+ assert(m_class.size()==m_threshold.size());
+ if(windowBuffer.size()){
+ outBuffer[x/down]=inBuffer[dimY/2][x];//initialize with original value (in case thresholds not met)
+ for(int iclass=0;iclass<m_class.size();++iclass){
+ if(100.0*(occurrence[m_class[iclass]])/windowBuffer.size()>m_threshold[iclass])
+ outBuffer[x/down]=m_class[iclass];
+ }
+ }
+ else
+ outBuffer[x/down]=m_noValue;
+ break;
+ }
+ case(MIXED):{
+ enum Type { BF=11, CF=12, MF=13, NF=20, W=30 };
+ double nBF=occurrence[BF];
+ double nCF=occurrence[CF];
+ double nMF=occurrence[MF];
+ double nNF=occurrence[NF];
+ double nW=occurrence[W];
+ if(windowBuffer.size()){
+ if((nBF+nCF+nMF)&&(nBF+nCF+nMF>=nNF+nW)){//forest
+ if(nBF/(nBF+nCF)>=0.75)
+ outBuffer[x/down]=BF;
+ else if(nCF/(nBF+nCF)>=0.75)
+ outBuffer[x/down]=CF;
+ else
+ outBuffer[x/down]=MF;
+ }
+ else{//non-forest
+ if(nW&&(nW>=nNF))
+ outBuffer[x/down]=W;
+ else
+ outBuffer[x/down]=NF;
+ }
+ }
+ else
+ outBuffer[x/down]=inBuffer[indexJ][indexI];
+ break;
+ }
+ default:
+ break;
}
}
+ progress=(1.0+y/down);
+ progress+=(outputVector.size());
+ progress/=outputVector.size();
+ pfnProgress(progress,pszMessage,pProgressArg);
//copy outBuffer to outputVector
- outputVector[y]=outBuffer;
+ outputVector[y/down]=outBuffer;
}
}
diff --git a/src/apps/pkclassify_svm.cc b/src/apps/pkclassify_svm.cc
index 3c86021..94139de 100644
--- a/src/apps/pkclassify_svm.cc
+++ b/src/apps/pkclassify_svm.cc
@@ -17,6 +17,7 @@ GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with pktools. If not, see <http://www.gnu.org/licenses/>.
***********************************************************************/
+
#include <vector>
#include <map>
#include <algorithm>
diff --git a/src/apps/pkcrop.cc b/src/apps/pkcrop.cc
index f17dc28..263f0ca 100644
--- a/src/apps/pkcrop.cc
+++ b/src/apps/pkcrop.cc
@@ -60,7 +60,7 @@ int main(int argc, char *argv[])
Optionpk<double> offset_opt("off", "offset", "output=scale*input+offset", 0);
Optionpk<string> otype_opt("ot", "otype", "Data type for output image ({Byte/Int16/UInt16/UInt32/Int32/Float32/Float64/CInt16/CInt32/CFloat32/CFloat64}). Empty string: inherit type from input image","");
Optionpk<string> oformat_opt("of", "oformat", "Output image format (see also gdal_translate). Empty string: inherit from input image");
- Optionpk<string> colorTable_opt("ct", "ct", "color table (file with 5 columns: id R G B ALFA (0: transparent, 255: solid)", "");
+ Optionpk<string> colorTable_opt("ct", "ct", "color table (file with 5 columns: id R G B ALFA (0: transparent, 255: solid)");
Optionpk<string> option_opt("co", "co", "options: NAME=VALUE [-co COMPRESS=LZW] [-co INTERLEAVE=BAND]");
Optionpk<short> flag_opt("f", "flag", "Flag value to put in image if out of bounds.", 0);
Optionpk<string> resample_opt("r", "resampling-method", "Resampling method (near: nearest neighbour, bilinear: bi-linear interpolation).", "near");
@@ -327,7 +327,7 @@ int main(int argc, char *argv[])
}
else if(imgReader.isGeoRef())
imgWriter.setProjection(imgReader.getProjection());
- if(colorTable_opt[0]!=""){
+ if(colorTable_opt.size()){
if(verbose_opt[0])
cout << "set colortable " << colorTable_opt[0] << endl;
assert(imgWriter.getDataType()==GDT_Byte);
diff --git a/src/apps/pkfilter.cc b/src/apps/pkfilter.cc
index f77fb21..7beadb0 100644
--- a/src/apps/pkfilter.cc
+++ b/src/apps/pkfilter.cc
@@ -50,8 +50,8 @@ int main(int argc,char **argv) {
Optionpk<bool> license_opt("lic","license","show license information",false);
Optionpk<bool> help_opt("h","help","shows this help info",false);
Optionpk<bool> todo_opt("\0","todo","",false);
- Optionpk<std::string> input_opt("i","input","input image file","");
- Optionpk<std::string> output_opt("o", "output", "Output image file", "");
+ Optionpk<std::string> input_opt("i","input","input image file");
+ Optionpk<std::string> output_opt("o", "output", "Output image file");
Optionpk<bool> disc_opt("c", "circular", "circular disc kernel for dilation and erosion", false);
Optionpk<double> angle_opt("a", "angle", "angle used for directional filtering in dilation.");
Optionpk<int> function_opt("f", "filter", "filter function (0: median, 1: variance, 2: min, 3: max, 4: sum, 5: mean, 6: minmax, 7: dilation, 8: erosion, 9: closing, 10: opening, 11: spatially homogeneous (central pixel must be identical to all other pixels within window), 12: SobelX edge detection in X, 13: SobelY edge detection in Y, 14: SobelXY, -14: SobelYX, 15: smooth, 16: density, 17: majority voting (only for classes), 18: forest aggregation (mixed), 19: smooth no data (mask) val [...]
@@ -60,12 +60,12 @@ int main(int argc,char **argv) {
Optionpk<int> dimZ_opt("dz", "dz", "filter kernel size in z (band or spectral dimension), must be odd (example: 3).. Set dz>0 if 1-D filter must be used in band domain");
Optionpk<short> class_opt("class", "class", "class value(s) to use for density, erosion, dilation, openening and closing, thresholding");
Optionpk<double> threshold_opt("t", "threshold", "threshold value(s) to use for threshold filter (one for each class)", 0);
- Optionpk<short> mask_opt("\0", "mask", "mask value(s) ");
+ Optionpk<short> mask_opt("m", "mask", "mask value(s) ");
Optionpk<std::string> tap_opt("tap", "tap", "text file containing taps used for spatial filtering (from ul to lr). Use dimX and dimY to specify tap dimensions in x and y. Leave empty for not using taps");
Optionpk<double> tapz_opt("tapz", "tapz", "taps used for spectral filtering");
- Optionpk<std::string> otype_opt("ot", "otype", "Data type for output image ({Byte/Int16/UInt16/UInt32/Int32/Float32/Float64/CInt16/CInt32/CFloat32/CFloat64}). Empty string: inherit type from input image", "");
+ Optionpk<std::string> otype_opt("ot", "otype", "Data type for output image ({Byte/Int16/UInt16/UInt32/Int32/Float32/Float64/CInt16/CInt32/CFloat32/CFloat64}). Empty string: inherit type from input image","");
Optionpk<string> oformat_opt("of", "oformat", "Output image format (see also gdal_translate). Empty string: inherit from input image");
- Optionpk<string> colorTable_opt("ct", "ct", "color table (file with 5 columns: id R G B ALFA (0: transparent, 255: solid)", "");
+ Optionpk<string> colorTable_opt("ct", "ct", "color table (file with 5 columns: id R G B ALFA (0: transparent, 255: solid)");
Optionpk<std::string> option_opt("co", "co", "options: NAME=VALUE [-co COMPRESS=LZW] [-co INTERLEAVE=BAND]");
Optionpk<short> down_opt("d", "down", "down sampling factor. Use value 1 for no downsampling)", 1);
Optionpk<short> verbose_opt("v", "verbose", "verbose mode if > 0", 0);
@@ -111,6 +111,7 @@ int main(int argc,char **argv) {
ImgReaderGdal input;
ImgWriterGdal output;
+ assert(input_opt.size());
input.open(input_opt[0]);
// output.open(output_opt[0],input);
GDALDataType theType=GDT_Unknown;
@@ -140,6 +141,7 @@ int main(int argc,char **argv) {
option_opt.push_back(theInterleave);
}
try{
+ assert(output_opt.size());
output.open(output_opt[0],(input.nrOfCol()+down_opt[0]-1)/down_opt[0],(input.nrOfRow()+down_opt[0]-1)/down_opt[0],input.nrOfBand(),theType,imageType,option_opt);
}
catch(string errorstring){
@@ -152,7 +154,13 @@ int main(int argc,char **argv) {
input.getGeoTransform(ulx,uly,deltaX,deltaY,rot1,rot2);
output.setGeoTransform(ulx,uly,deltaX*down_opt[0],deltaY*down_opt[0],rot1,rot2);
}
- if(input.getColorTable()!=NULL)
+ if(colorTable_opt.size()){
+ if(verbose_opt[0])
+ cout << "set colortable " << colorTable_opt[0] << endl;
+ assert(output.getDataType()==GDT_Byte);
+ output.setColorTable(colorTable_opt[0]);
+ }
+ else if(input.getColorTable()!=NULL)
output.setColorTable(input.getColorTable());
// if(input.isGeoRef()){
@@ -162,7 +170,7 @@ int main(int argc,char **argv) {
Filter2d::Filter2d filter2d;
Filter filter1d;
- if(verbose_opt[0])
+ if(class_opt.size()&&verbose_opt[0])
std::cout<< "class values: ";
for(int iclass=0;iclass<class_opt.size();++iclass){
if(!dimZ_opt.size())
@@ -211,7 +219,7 @@ int main(int argc,char **argv) {
filter1d.doit(input,output,down_opt[0]);
}
else{
- if(colorTable_opt[0]!="")
+ if(colorTable_opt.size())
output.setColorTable(colorTable_opt[0]);
switch(function_opt[0]){
case(Filter2d::MEDIAN):
diff --git a/src/apps/pklas2img.cc b/src/apps/pklas2img.cc
index fe26312..4233052 100644
--- a/src/apps/pklas2img.cc
+++ b/src/apps/pklas2img.cc
@@ -50,9 +50,10 @@ int main(int argc,char **argv) {
Optionpk<double> hThreshold_opt("ht", "maxHeight", "initial and maximum height threshold for progressive morphological filtering (e.g., -ht 0.2 -ht 2.5)", 0.2);
Optionpk<short> maxIter_opt("\0", "maxIter", "Maximum number of iterations in post filter (default is 100)", 100.0);
Optionpk<short> nbin_opt("nb", "nbin", "Number of percentile bins for calculating profile (=number of output bands) (default is 10)", 10.0);
- Optionpk<unsigned short> returns_opt("r", "returns", "number(s) of returns to include (use -r -1 for last return only). Default is 0 (include all returns)", 0);
- Optionpk<string> composite_opt("c", "composite", "composite for multiple points in cell (min, max, median, mean, sum, first, last, profile, number (point density)). Default is last (overwrite cells with latest point", "last");
- Optionpk<string> filter_opt("fir", "filter", "filter las points (single,multiple,all). Default is all", "all");
+ Optionpk<unsigned short> returns_opt("r", "returns", "number(s) of returns to include");
+ Optionpk<unsigned short> classes_opt("c", "classes", "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, number (point density)). Default is last (overwrite cells with latest point", "last");
+ Optionpk<string> filter_opt("fir", "filter", "filter las points (last,single,multiple,all). Default is all", "all");
Optionpk<string> postFilter_opt("pf", "pfilter", "post processing filter (etew_min,promorph (progressive morphological filter),bunting (adapted promorph),open,close,none) . Default is none", "none");
Optionpk<short> dimx_opt("\0", "dimX", "Dimension X of postFilter (default is 3)", 3);
Optionpk<short> dimy_opt("\0", "dimY", "Dimension Y of postFilter (default is 3)", 3);
@@ -94,6 +95,7 @@ int main(int argc,char **argv) {
maxIter_opt.retrieveOption(argc,argv);
nbin_opt.retrieveOption(argc,argv);
returns_opt.retrieveOption(argc,argv);
+ classes_opt.retrieveOption(argc,argv);
composite_opt.retrieveOption(argc,argv);
filter_opt.retrieveOption(argc,argv);
postFilter_opt.retrieveOption(argc,argv);
@@ -256,8 +258,10 @@ int main(int argc,char **argv) {
//set bounding filter
// lasReader.addBoundsFilter(minULX,maxULY,maxLRX,minLRY);
//set returns filter
- if(returns_opt[0])
+ if(returns_opt.size())
lasReader.addReturnsFilter(returns_opt);
+ if(classes_opt.size())
+ lasReader.addClassFilter(classes_opt);
lasReader.setFilters();
if(attribute_opt[0]!="z"){
@@ -295,6 +299,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;
+ }
double dcol,drow;
outputWriter.geo2image(thePoint.GetX(),thePoint.GetY(),dcol,drow);
int icol=static_cast<int>(dcol);
@@ -454,12 +462,10 @@ int main(int argc,char **argv) {
try{
theFilter.morphology(outputData,currentOutput,Filter2d::ERODE,dimx,dimy,disc_opt[0],maxSlope_opt[0]);
theFilter.morphology(currentOutput,outputData,Filter2d::DILATE,dimx,dimy,disc_opt[0],maxSlope_opt[0]);
- // if(postFilter_opt[0]=="bunting"){//todo: implement doit in Filter2d on Vector2d
- // theFilter.doit(outputData,currentOutput,Filter2d::MEDIAN,dimx,dimy,1,disc_opt[0]);
- // filter2d.doit(input,output,Filter2d::MEDIAN,dimX_opt[0],dimY_opt[0],down_opt[0],disc_opt[0]);
-
- // outputData=currentOutput;
- // }
+ if(postFilter_opt[0]=="bunting"){//todo: implement doit in Filter2d on Vector2d
+ theFilter.doit(outputData,currentOutput,Filter2d::MEDIAN,dimx,dimy,1,disc_opt[0]);
+ outputData=currentOutput;
+ }
}
catch(std::string errorString){
cout << errorString << endl;
diff --git a/src/fileclasses/FileReaderLas.cc b/src/fileclasses/FileReaderLas.cc
index 25f1d42..aceaac1 100644
--- a/src/fileclasses/FileReaderLas.cc
+++ b/src/fileclasses/FileReaderLas.cc
@@ -22,6 +22,30 @@ along with pktools. If not, see <http://www.gnu.org/licenses/>.
#include <fstream>
#include "FileReaderLas.h"
//---------------------------------------------------------------------------
+LastReturnFilter::LastReturnFilter( ) : liblas::FilterI(eInclusion) {}
+
+bool LastReturnFilter::filter(const liblas::Point& p)
+{
+
+ // If the GetReturnNumber equals the GetNumberOfReturns,
+ // we're a last return
+
+ bool output = false;
+ if (p.GetReturnNumber() == p.GetNumberOfReturns())
+ {
+ output = true;
+ }
+
+ // If the type is switched to eExclusion, we'll throw out all last returns.
+ if (GetType() == eExclusion && output == true)
+ {
+ output = false;
+ } else {
+ output = true;
+ }
+ return output;
+}
+
FileReaderLas::FileReaderLas(void)
{
m_reader=NULL;
@@ -138,3 +162,18 @@ void FileReaderLas::addReturnsFilter(std::vector<unsigned short> const& returns)
m_filters.push_back(liblas::FilterPtr(return_filter));
}
+void FileReaderLas::addClassFilter(std::vector<unsigned short> const& classes){
+
+ std::vector<liblas::FilterPtr> filters;
+ std::vector<liblas::Classification> theClasses;
+ for(int iclass=0;iclass<classes.size();++iclass){
+ liblas::Classification aClass(classes[iclass]);
+ theClasses.push_back(aClass);
+ }
+ liblas::FilterPtr class_filter = liblas::FilterPtr(new liblas::ClassificationFilter(theClasses));
+ // eInclusion means to keep the classes that match. eExclusion would
+ // throw out those that matched
+ class_filter->SetType(liblas::FilterI::eInclusion);
+ m_filters.push_back(class_filter);
+}
+
diff --git a/src/fileclasses/FileReaderLas.h b/src/fileclasses/FileReaderLas.h
index 93023dd..8435d01 100644
--- a/src/fileclasses/FileReaderLas.h
+++ b/src/fileclasses/FileReaderLas.h
@@ -25,6 +25,17 @@ along with pktools. If not, see <http://www.gnu.org/licenses/>.
#include "liblas/liblas.hpp"
//--------------------------------------------------------------------------
+class LastReturnFilter: public liblas::FilterI
+{
+public:
+ LastReturnFilter();
+ bool filter(const liblas::Point& point);
+
+private:
+ LastReturnFilter(LastReturnFilter const& other);
+ LastReturnFilter& operator=(LastReturnFilter const& rhs);
+};
+
class FileReaderLas
{
public:
@@ -47,6 +58,7 @@ public:
liblas::Point const& readPointAt(std::size_t n){m_reader->ReadPointAt(n);return m_reader->GetPoint();};
// void addBoundsFilter(double ulx, double uly, double lrx, double lry);
void addReturnsFilter(std::vector<unsigned short> const& returns);
+ void addClassFilter(std::vector<unsigned short> const& classes);
void setFilters(const std::vector<liblas::FilterPtr>& filters){m_filters=filters;setFilters();};
void setFilters(){m_reader->SetFilters(m_filters);};
protected:
--
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