[pktools] 323/375: working on filter in spectral domain
Bas Couwenberg
sebastic at xs4all.nl
Wed Dec 3 21:54:26 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 872e3f286852671602007807f9f536a470157033
Author: Pieter Kempeneers <kempenep at gmail.com>
Date: Sun Sep 21 11:05:51 2014 -0700
working on filter in spectral domain
---
src/algorithms/Filter.cc | 87 +++++++++++++++++++++
src/algorithms/Filter.h | 199 +++++++++++++++++++++++++++++++++++++++++++++++
src/apps/pkfilter.cc | 10 ++-
3 files changed, 294 insertions(+), 2 deletions(-)
diff --git a/src/algorithms/Filter.cc b/src/algorithms/Filter.cc
index eb2ddb4..58080bd 100644
--- a/src/algorithms/Filter.cc
+++ b/src/algorithms/Filter.cc
@@ -326,6 +326,93 @@ void filter::Filter::filter(const ImgReaderGdal& input, ImgWriterGdal& output, s
}
}
+void filter::Filter::stat(const ImgReaderGdal& input, ImgWriterGdal& output, const std::string& method, short down, int offset)
+{
+ Vector2d<double> lineInput(input.nrOfBand(),input.nrOfCol());
+ assert(output.nrOfCol()==(input.nrOfCol()+down-1)/down);
+ vector<double> lineOutput(output.nrOfCol());
+ statfactory::StatFactory stat;
+ const char* pszMessage;
+ void* pProgressArg=NULL;
+ GDALProgressFunc pfnProgress=GDALTermProgress;
+ double progress=0;
+ pfnProgress(progress,pszMessage,pProgressArg);
+ for(int y=0;y<input.nrOfRow();++y){
+ for(int iband=0;iband<input.nrOfBand();++iband)
+ input.readData(lineInput[iband],GDT_Float64,y,iband);
+ vector<double> pixelInput(input.nrOfBand());
+ for(int x=0;x<input.nrOfCol();++x){
+ pixelInput=lineInput.selectCol(x);
+ switch(getFilterType(method)){
+ case(filter::median):
+ lineOutput[(x-offset+down-1)/down]=stat.median(pixelInput);
+ break;
+ case(filter::min):
+ lineOutput[(x-offset+down-1)/down]=stat.mymin(pixelInput);
+ break;
+ case(filter::max):
+ lineOutput[(x-offset+down-1)/down]=stat.mymax(pixelInput);
+ break;
+ case(filter::sum):
+ lineOutput[(x-offset+down-1)/down]=sqrt(stat.sum(pixelInput));
+ break;
+ case(filter::var):
+ lineOutput[(x-offset+down-1)/down]=stat.var(pixelInput);
+ break;
+ case(filter::mean):
+ lineOutput[(x-offset+down-1)/down]=stat.mean(pixelInput);
+ break;
+ default:
+ std::string errorString="method not supported";
+ throw(errorString);
+ break;
+ }
+ }
+ try{
+ output.writeData(lineOutput,GDT_Float64,y);
+ }
+ catch(string errorstring){
+ cerr << errorstring << "in line " << y << endl;
+ }
+ progress=(1.0+y)/output.nrOfRow();
+ pfnProgress(progress,pszMessage,pProgressArg);
+ }
+}
+
+void filter::Filter::filter(const ImgReaderGdal& input, ImgWriterGdal& output, const std::string& method, int dim, short down, int offset)
+{
+ Vector2d<double> lineInput(input.nrOfBand(),input.nrOfCol());
+ Vector2d<double> lineOutput;
+ const char* pszMessage;
+ void* pProgressArg=NULL;
+ GDALProgressFunc pfnProgress=GDALTermProgress;
+ double progress=0;
+ pfnProgress(progress,pszMessage,pProgressArg);
+ for(int y=0;y<input.nrOfRow();++y){
+ for(int iband=0;iband<input.nrOfBand();++iband)
+ input.readData(lineInput[iband],GDT_Float64,y,iband);
+ vector<double> pixelInput(input.nrOfBand());
+ vector<double> pixelOutput;
+ for(int x=0;x<input.nrOfCol();++x){
+ pixelInput=lineInput.selectCol(x);
+ filter(pixelInput,pixelOutput,method,dim,down,offset);
+ lineOutput.resize(pixelOutput.size(),input.nrOfCol());
+ for(int iband=0;iband<pixelOutput.size();++iband)
+ lineOutput[iband][x]=pixelOutput[iband];
+ }
+ for(int iband=0;iband<input.nrOfBand();++iband){
+ try{
+ output.writeData(lineOutput[iband],GDT_Float64,y,iband);
+ }
+ catch(string errorstring){
+ cerr << errorstring << "in band " << iband << ", line " << y << endl;
+ }
+ }
+ progress=(1.0+y)/output.nrOfRow();
+ pfnProgress(progress,pszMessage,pProgressArg);
+ }
+}
+
double filter::Filter::getCentreWavelength(const std::vector<double> &wavelengthIn, const Vector2d<double>& srf, const std::string& interpolationType, double delta, bool verbose)
{
assert(srf.size()==2);//[0]: wavelength, [1]: response function
diff --git a/src/algorithms/Filter.h b/src/algorithms/Filter.h
index dc4cc12..b8172e8 100644
--- a/src/algorithms/Filter.h
+++ b/src/algorithms/Filter.h
@@ -58,12 +58,15 @@ public:
void pushClass(short theClass=1){m_class.push_back(theClass);};
void pushMask(short theMask=0){m_mask.push_back(theMask);};
template<class T> void filter(const std::vector<T>& input, std::vector<T>& output, int down=1, int offset=0);
+ template<class T> void filter(const std::vector<T>& input, std::vector<T>& output, const std::string& method, int dim, int down=1, int offset=0);
template<class T> void smooth(const std::vector<T>& input, std::vector<T>& output, short dim, int down=1, int offset=0);
template<class T> void filter(T* input, int inputSize, std::vector<T>& output, int down=1, int offset=0);
template<class T> void smooth(T* input, int inputSize, std::vector<T>& output, short dim, int down=1, int offset=0);
template<class T> void morphology(const std::vector<T>& input, std::vector<T>& output, const std::string& method, int dim, short down=1, int offset=0, bool verbose=false);
void morphology(const ImgReaderGdal& input, ImgWriterGdal& output, const std::string& method, int dim, short down=1, int offset=0, short verbose=0);
void filter(const ImgReaderGdal& input, ImgWriterGdal& output, short down=1, int offset=0);
+ void stat(const ImgReaderGdal& input, ImgWriterGdal& output, const std::string& method, short down=1, int offset=0);
+ void filter(const ImgReaderGdal& input, ImgWriterGdal& output, const std::string& method, int dim, short down=1, int offset=0);
void smooth(const ImgReaderGdal& input, ImgWriterGdal& output, short dim, short down=1, int offset=0);
double getCentreWavelength(const std::vector<double> &wavelengthIn, const Vector2d<double>& srf, const std::string& interpolationType, double delta=1.0, bool verbose=false);
template<class T> double applySrf(const std::vector<double> &wavelengthIn, const std::vector<T>& input, const Vector2d<double>& srf, const std::string& interpolationType, T& output, double delta=1.0, bool normalize=false, bool verbose=false);
@@ -447,6 +450,202 @@ template<class T> void Filter::filter(const std::vector<T>& input, std::vector<T
}
}
+ template<class T> void Filter::filter(const std::vector<T>& input, std::vector<T>& output, const std::string& method, int dim, int down, int offset)
+{
+ bool verbose=false;
+ assert(dim);
+ output.resize((input.size()-offset+down-1)/down);
+ int i=0;
+ statfactory::StatFactory stat;
+ std::vector<T> statBuffer;
+ short binValue=0;
+ //start: extend input with mirrored version of itself
+ for(i=offset;i<dim/2;++i){
+ binValue=0;
+ for(int iclass=0;iclass<m_class.size();++iclass){
+ if(input[i]==m_class[iclass]){
+ binValue=m_class[0];
+ break;
+ }
+ }
+ if(m_class.size())
+ statBuffer.push_back(binValue);
+ else
+ statBuffer.push_back(input[i]);
+ for(int t=1;t<=dim/2;++t){
+ binValue=0;
+ for(int iclass=0;iclass<m_class.size();++iclass){
+ if(input[i+t]==m_class[iclass]){
+ binValue=m_class[0];
+ break;
+ }
+ }
+ if(m_class.size()){
+ statBuffer.push_back(binValue);
+ statBuffer.push_back(binValue);
+ }
+ else{
+ statBuffer.push_back(input[i+t]);
+ statBuffer.push_back(input[i+t]);
+ }
+ }
+ assert(statBuffer.size()==dim);
+ if((i-offset)%down){
+ statBuffer.clear();
+ continue;
+ }
+ switch(getFilterType(method)){
+ case(filter::median):
+ output[(i-offset+down-1)/down]=stat.median(statBuffer);
+ break;
+ case(filter::min):
+ output[(i-offset+down-1)/down]=stat.mymin(statBuffer);
+ break;
+ case(filter::max):
+ output[(i-offset+down-1)/down]=stat.mymax(statBuffer);
+ break;
+ case(filter::sum):
+ output[(i-offset+down-1)/down]=sqrt(stat.sum(statBuffer));
+ break;
+ case(filter::var):
+ output[(i-offset+down-1)/down]=stat.var(statBuffer);
+ break;
+ case(filter::mean):
+ output[(i-offset+down-1)/down]=stat.mean(statBuffer);
+ break;
+ default:
+ std::string errorString="method not supported";
+ throw(errorString);
+ break;
+ }
+ if(verbose){
+ std::cout << "buffer: ";
+ for(int ibuf=0;ibuf<statBuffer.size();++ibuf)
+ std::cout << statBuffer[ibuf] << " ";
+ std::cout << "->" << output[(i-offset+down-1)/down] << std::endl;
+ }
+ }
+ //main
+ statBuffer.clear();
+ for(i=offset+dim/2;i<input.size()-dim/2;++i){
+ binValue=0;
+ for(int t=0;t<dim;++t){
+ for(int iclass=0;iclass<m_class.size();++iclass){
+ if(input[i-dim/2+t]==m_class[iclass]){
+ binValue=m_class[0];
+ break;
+ }
+ }
+ if(m_class.size())
+ statBuffer.push_back(binValue);
+ else
+ statBuffer.push_back(input[i-dim/2+t]);
+ }
+ assert(statBuffer.size()==dim);
+ if((i-offset)%down){
+ statBuffer.clear();
+ continue;
+ }
+ switch(getFilterType(method)){
+ case(filter::median):
+ output[(i-offset+down-1)/down]=stat.median(statBuffer);
+ break;
+ case(filter::min):
+ output[(i-offset+down-1)/down]=stat.mymin(statBuffer);
+ break;
+ case(filter::max):
+ output[(i-offset+down-1)/down]=stat.mymax(statBuffer);
+ break;
+ case(filter::sum):
+ output[(i-offset+down-1)/down]=sqrt(stat.sum(statBuffer));
+ break;
+ case(filter::var):
+ output[(i-offset+down-1)/down]=stat.var(statBuffer);
+ break;
+ case(filter::mean):
+ output[(i-offset+down-1)/down]=stat.mean(statBuffer);
+ break;
+ default:
+ std::string errorString="method not supported";
+ throw(errorString);
+ break;
+ }
+ if(verbose){
+ std::cout << "buffer: ";
+ for(int ibuf=0;ibuf<statBuffer.size();++ibuf)
+ std::cout << statBuffer[ibuf] << " ";
+ std::cout << "->" << output[(i-offset+down-1)/down] << std::endl;
+ }
+ statBuffer.clear();
+ }
+ //end: extend input with mirrored version of itself
+ for(i=input.size()-dim/2;i<input.size();++i){
+ binValue=0;
+ for(int iclass=0;iclass<m_class.size();++iclass){
+ if(input[i]==m_class[iclass]){
+ binValue=m_class[0];
+ break;
+ }
+ }
+ if(m_class.size())
+ statBuffer.push_back(binValue);
+ else
+ statBuffer.push_back(input[i]);
+ for(int t=1;t<=dim/2;++t){
+ binValue=0;
+ for(int iclass=0;iclass<m_class.size();++iclass){
+ if(input[i-t]==m_class[iclass]){
+ binValue=m_class[0];
+ break;
+ }
+ }
+ if(m_class.size()){
+ statBuffer.push_back(binValue);
+ statBuffer.push_back(binValue);
+ }
+ else{
+ statBuffer.push_back(input[i-t]);
+ statBuffer.push_back(input[i-t]);
+ }
+ }
+ if((i-offset)%down){
+ statBuffer.clear();
+ continue;
+ }
+ switch(getFilterType(method)){
+ case(filter::median):
+ output[(i-offset+down-1)/down]=stat.median(statBuffer);
+ break;
+ case(filter::min):
+ output[(i-offset+down-1)/down]=stat.mymin(statBuffer);
+ break;
+ case(filter::max):
+ output[(i-offset+down-1)/down]=stat.mymax(statBuffer);
+ break;
+ case(filter::sum):
+ output[(i-offset+down-1)/down]=sqrt(stat.sum(statBuffer));
+ break;
+ case(filter::var):
+ output[(i-offset+down-1)/down]=stat.var(statBuffer);
+ break;
+ case(filter::mean):
+ output[(i-offset+down-1)/down]=stat.mean(statBuffer);
+ break;
+ default:
+ std::string errorString="method not supported";
+ throw(errorString);
+ break;
+ }
+ if(verbose){
+ std::cout << "buffer: ";
+ for(int ibuf=0;ibuf<statBuffer.size();++ibuf)
+ std::cout << statBuffer[ibuf] << " ";
+ std::cout << "->" << output[(i-offset+down-1)/down] << std::endl;
+ }
+ }
+}
+
+ //todo: this function is redundant, can be incorporated within filter
template<class T> void Filter::morphology(const std::vector<T>& input, std::vector<T>& output, const std::string& method, int dim, short down, int offset, bool verbose)
{
assert(dim);
diff --git a/src/apps/pkfilter.cc b/src/apps/pkfilter.cc
index 11bc602..d8f04cd 100644
--- a/src/apps/pkfilter.cc
+++ b/src/apps/pkfilter.cc
@@ -188,8 +188,10 @@ int main(int argc,char **argv) {
int nband=fwhm_opt.size()? fwhm_opt.size():srf_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],nband,theType,imageType,option_opt);
}
- else
- 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);
+ else{
+ int nband=(dimZ_opt[0]==1)? 1 : input.nrOfBand();
+ output.open(output_opt[0],(input.nrOfCol()+down_opt[0]-1)/down_opt[0],(input.nrOfRow()+down_opt[0]-1)/down_opt[0],nband,theType,imageType,option_opt);
+ }
}
catch(string errorstring){
cout << errorstring << endl;
@@ -707,6 +709,10 @@ int main(int argc,char **argv) {
std::cout << "classes set" << std::endl;
default:
if(dimZ_opt.size()){
+ if(dimZ_opt.size()==1)
+ filter1d.stat(input,output,method_opt[0],down_opt[0]);
+ else
+ filter1d.filter(input,output,method_opt[0],dimZ_opt[0],down_opt[0]);
}
else
filter2d.doit(input,output,method_opt[0],dimX_opt[0],dimY_opt[0],down_opt[0],disc_opt[0]);
--
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