[pktools] 171/375: support for non shape files extended, problem in pkdiff when output is defined remains...

Bas Couwenberg sebastic at xs4all.nl
Wed Dec 3 21:54:11 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 e25b014405ca1c5bdb35cc3e50318910a9ccc9fb
Author: Pieter Kempeneers <kempenep at gmail.com>
Date:   Sat Dec 28 22:47:28 2013 +0100

    support for non shape files extended, problem in pkdiff when output is defined remains...
---
 src/apps/pkclassify_nn.cc        |  16 +-
 src/apps/pkclassify_svm.cc       |  15 +-
 src/apps/pkdiff.cc               | 665 ++++++++++++++++++++-------------------
 src/apps/pkeditogr.cc            |  21 +-
 src/apps/pkreclass.cc            |  13 +-
 src/imageclasses/ImgWriterOgr.cc |   2 +-
 6 files changed, 397 insertions(+), 335 deletions(-)

diff --git a/src/apps/pkclassify_nn.cc b/src/apps/pkclassify_nn.cc
index 147d259..e369096 100644
--- a/src/apps/pkclassify_nn.cc
+++ b/src/apps/pkclassify_nn.cc
@@ -60,7 +60,7 @@ int main(int argc, char *argv[])
   Optionpk<unsigned short> bag_opt("bag", "bag", "Number of bootstrap aggregations (default is no bagging: 1)", 1);
   Optionpk<int> bagSize_opt("bs", "bsize", "Percentage of features used from available training features for each bootstrap aggregation (one size for all classes, or a different size for each class respectively", 100);
   Optionpk<string> classBag_opt("cb", "classbag", "output for each individual bootstrap aggregation (default is blank)"); 
-  Optionpk<string> mask_opt("m", "mask", "mask image (see also mvalue option (default is no mask)"); 
+  Optionpk<string> mask_opt("m", "mask", "mask image (support for single mask only, see also msknodata option)"); 
   Optionpk<short> msknodata_opt("msknodata", "msknodata", "mask value(s) not to consider for classification (use negative values if only these values should be taken into account). Values will be taken over in classification image. Default is 0", 0);
   Optionpk<unsigned short> nodata_opt("nodata", "nodata", "nodata value to put where image is masked as nodata", 0);
   Optionpk<string> output_opt("o", "output", "output classification image"); 
@@ -573,7 +573,17 @@ int main(int argc, char *argv[])
     GDALProgressFunc pfnProgress=GDALTermProgress;
     float progress=0;
   //-------------------------------- open image file ------------------------------------
-  if(input_opt[0].find(".shp")==string::npos){
+  bool refIsRaster=false;
+  ImgReaderOgr imgReaderOgr;
+  try{
+    imgReaderOgr.open(input_opt[0]);
+    imgReaderOgr.close();
+  }
+  catch(string errorString){
+    refIsRaster=true;
+  }
+  if(refIsRaster){
+  // if(input_opt[0].find(".shp")==string::npos){
     ImgReaderGdal testImage;
     try{
       if(verbose_opt[0]>=1)
@@ -957,7 +967,7 @@ int main(int argc, char *argv[])
         assert(output_opt.size()==input_opt.size());
       if(verbose_opt[0])
         cout << "opening img reader " << input_opt[ivalidation] << endl;
-      ImgReaderOgr imgReaderOgr(input_opt[ivalidation]);
+      imgReaderOgr.open(input_opt[ivalidation]);
       ImgWriterOgr imgWriterOgr;
 
       if(output_opt.size()){
diff --git a/src/apps/pkclassify_svm.cc b/src/apps/pkclassify_svm.cc
index 96ee9e1..8a1512c 100644
--- a/src/apps/pkclassify_svm.cc
+++ b/src/apps/pkclassify_svm.cc
@@ -77,7 +77,7 @@ int main(int argc, char *argv[])
   Optionpk<unsigned short> bag_opt("bag", "bag", "Number of bootstrap aggregations", 1);
   Optionpk<int> bagSize_opt("bs", "bsize", "Percentage of features used from available training features for each bootstrap aggregation (one size for all classes, or a different size for each class respectively", 100);
   Optionpk<string> classBag_opt("cb", "classbag", "output for each individual bootstrap aggregation");
-  Optionpk<string> mask_opt("m", "mask", "mask image (see also mvalue option"); 
+  Optionpk<string> mask_opt("m", "mask", "mask image (support for single mask only, see also msknodata option)"); 
   Optionpk<short> msknodata_opt("msknodata", "msknodata", "mask value(s) not to consider for classification (use negative values if only these values should be taken into account). Values will be taken over in classification image.", 0);
   Optionpk<unsigned short> nodata_opt("nodata", "nodata", "nodata value to put where image is masked as nodata", 0);
   Optionpk<string> output_opt("o", "output", "output classification image"); 
@@ -561,7 +561,16 @@ int main(int argc, char *argv[])
   if(!verbose_opt[0])
     pfnProgress(progress,pszMessage,pProgressArg);
   //-------------------------------- open image file ------------------------------------
-  if(input_opt[0].find(".shp")==string::npos){
+  bool refIsRaster=false;
+  ImgReaderOgr imgReaderOgr;
+  try{
+    imgReaderOgr.open(input_opt[0]);
+    imgReaderOgr.close();
+  }
+  catch(string errorString){
+    refIsRaster=true;
+  }
+  if(refIsRaster){
     ImgReaderGdal testImage;
     try{
       if(verbose_opt[0]>=1)
@@ -956,7 +965,7 @@ int main(int argc, char *argv[])
 	assert(output_opt.size()==input_opt.size());
       if(verbose_opt[0])
         std::cout << "opening img reader " << input_opt[ivalidation] << std::endl;
-      ImgReaderOgr imgReaderOgr(input_opt[ivalidation]);
+      imgReaderOgr.open(input_opt[ivalidation]);
       ImgWriterOgr imgWriterOgr;
 
       if(output_opt.size()){
diff --git a/src/apps/pkdiff.cc b/src/apps/pkdiff.cc
index 23d4139..ea0a597 100644
--- a/src/apps/pkdiff.cc
+++ b/src/apps/pkdiff.cc
@@ -97,8 +97,15 @@ int main(int argc, char *argv[])
     cout << endl;
   }
 
-  assert(input_opt.size());
-  assert(reference_opt.size());
+  if(input_opt.empty()){
+    std::cerr << "No input file provided (use option -i). Use --help for help information" << std::endl;
+    exit(0);
+  }
+  if(reference_opt.empty()){
+    std::cerr << "No reference file provided (use option -ref). Use --help for help information" << std::endl;
+    exit(0);
+  }
+
   if(mask_opt.size())
     while(mask_opt.size()<input_opt.size())
       mask_opt.push_back(mask_opt[0]);
@@ -188,15 +195,25 @@ int main(int argc, char *argv[])
   double progress=0;
   if(!verbose_opt[0])
     pfnProgress(progress,pszMessage,pProgressArg);
-  if(reference_opt[0].find(".shp")!=string::npos){
+
+  bool refIsRaster=false;
+  ImgReaderOgr referenceReaderOgr;
+  try{
+    referenceReaderOgr.open(reference_opt[0]);
+    referenceReaderOgr.close();
+  }
+  catch(string errorString){
+    refIsRaster=true;
+  }
+  // if(reference_opt[0].find(".shp")!=string::npos){
+  if(!refIsRaster){
     for(int iinput=0;iinput<input_opt.size();++iinput){
       if(output_opt.size())
         assert(reference_opt.size()==output_opt.size());
       for(int iref=0;iref<reference_opt.size();++iref){
         if(verbose_opt[0])
           cout << "reference is " << reference_opt[iref] << endl;
-        assert(reference_opt[iref].find(".shp")!=string::npos);
-        ImgReaderOgr referenceReader;
+        // assert(reference_opt[iref].find(".shp")!=string::npos);
         try{
           inputReader.open(input_opt[iinput]);//,imagicX_opt[0],imagicY_opt[0]);
           if(mask_opt.size()){
@@ -204,7 +221,7 @@ int main(int argc, char *argv[])
             assert(inputReader.nrOfCol()==maskReader.nrOfCol());
             assert(inputReader.nrOfRow()==maskReader.nrOfRow());
           }
-          referenceReader.open(reference_opt[iref]);
+          referenceReaderOgr.open(reference_opt[iref]);
         }
         catch(string error){
           cerr << error << endl;
@@ -214,308 +231,324 @@ int main(int argc, char *argv[])
           referenceRange=inputRange;
 
         ImgWriterOgr ogrWriter;
-        OGRLayer *writeLayer;
-        if(output_opt.size()){
-          if(verbose_opt[0])
-            cout << "creating output vector file " << output_opt[0] << endl;
-          assert(output_opt[0].find(".shp")!=string::npos);
-          try{
-            ogrWriter.open(output_opt[iref],ogrformat_opt[0]);
-          }
-          catch(string error){
-            cerr << error << endl;
-            exit(1);
-          }
-          char     **papszOptions=NULL;
-          string layername=output_opt[0].substr(0,output_opt[0].find(".shp"));
-          if(verbose_opt[0])
-            cout << "creating layer: " << layername << endl;
-          if(ogrWriter.createLayer(layername, "EPSG:3035", wkbPoint, papszOptions)==NULL)
-            cout << "Error: create layer failed!" << endl;
-          else if(verbose_opt[0])
-            cout << "created layer" << endl;
-          if(verbose_opt[0])
-            cout << "copy fields from " << reference_opt[iref] << endl;
-          ogrWriter.copyFields(referenceReader);
-          //create extra field for classified label
-          short theDim=boundary_opt[0];
-          for(int windowJ=-theDim/2;windowJ<(theDim+1)/2;++windowJ){
-            for(int windowI=-theDim/2;windowI<(theDim+1)/2;++windowI){
-              if(disc_opt[0]&&(windowI*windowI+windowJ*windowJ>(theDim/2)*(theDim/2)))
-                continue;
-              ostringstream fs;
-              if(theDim>1)
-                fs << labelclass_opt[0] << "_" << windowJ << "_" << windowI;
-              else
-                fs << labelclass_opt[0];
-              if(verbose_opt[0])
-                cout << "creating field " << fs.str() << endl;
-              ogrWriter.createField(fs.str(),OFTInteger);
-            }
-          }
-          writeLayer=ogrWriter.getDataSource()->GetLayer(0);
-        }
-        OGRLayer  *readLayer;
-        readLayer = referenceReader.getDataSource()->GetLayer(0);
-        readLayer->ResetReading();
-        OGRFeature *readFeature;
-        int isample=0;
-        while( (readFeature = readLayer->GetNextFeature()) != NULL ){
-          if(verbose_opt[0])
-            cout << "sample " << ++isample << endl;
-          //get x and y from readFeature
-          double x,y;
-          OGRGeometry *poGeometry;
-          OGRPolygon readPolygon;
-          OGRPoint centroidPoint;
-          OGRPoint *poPoint;
-          poGeometry = readFeature->GetGeometryRef();
-          // assert( poGeometry != NULL && wkbFlatten(poGeometry->getGeometryType()) == wkbPoint );
-          if(poGeometry==NULL)
-            continue;
-          else if(wkbFlatten(poGeometry->getGeometryType()) == wkbPolygon){
-            readPolygon = *((OGRPolygon *) poGeometry);
-            readPolygon.Centroid(&centroidPoint);
-            poPoint=¢roidPoint;
-          }
-          else if(wkbFlatten(poGeometry->getGeometryType()) == wkbPoint )
-            poPoint = (OGRPoint *) poGeometry;
-          else{
-            std::cerr << "Warning: skipping feature (not of type point or polygon)" << std::endl;
-            continue;
-          }
-          x=poPoint->getX();
-          y=poPoint->getY();
-          short inputValue;
-          vector<short> inputValues;
-          bool isHomogeneous=true;
-          short maskValue;
-          short outputValue;
-          //read referenceValue from feature
-          short referenceValue;
-          string referenceClassName;
-          if(classValueMap.size()){
-            referenceClassName=readFeature->GetFieldAsString(readFeature->GetFieldIndex(labelref_opt[0].c_str()));
-            referenceValue=classValueMap[referenceClassName];
-          }
-          else
-            referenceValue=readFeature->GetFieldAsInteger(readFeature->GetFieldIndex(labelref_opt[0].c_str()));
-          if(verbose_opt[0])
-            cout << "reference value: " << referenceValue << endl;
-          bool pixelFlagged=false;
-          bool maskFlagged=false;
-          for(int iflag=0;iflag<nodata_opt.size();++iflag){
-            if(referenceValue==nodata_opt[iflag])
-              pixelFlagged=true;
-          }
-          if(pixelFlagged)
-            continue;
-          double i_centre,j_centre;
-          //input reader is georeferenced!
-          inputReader.geo2image(x,y,i_centre,j_centre);
-          //       else{
-          //         i_centre=x;
-          //         j_centre=y;
-          //       }
-          //nearest neighbour
-          j_centre=static_cast<int>(j_centre);
-          i_centre=static_cast<int>(i_centre);
-          //check if j_centre is out of bounds
-          if(static_cast<int>(j_centre)<0||static_cast<int>(j_centre)>=inputReader.nrOfRow())
-            continue;
-          //check if i_centre is out of bounds
-          if(static_cast<int>(i_centre)<0||static_cast<int>(i_centre)>=inputReader.nrOfCol())
-            continue;
-          OGRFeature *writeFeature;
-          if(output_opt.size()){
-            writeFeature = OGRFeature::CreateFeature(writeLayer->GetLayerDefn());
-            if(verbose_opt[0])
-              cout << "copying fields from " << reference_opt[0] << endl;
-            if(writeFeature->SetFrom(readFeature)!= OGRERR_NONE)
-              cerr << "writing feature failed" << endl;
-          }
-          bool windowAllFlagged=true;
-          bool windowHasFlag=false;
-          short theDim=boundary_opt[0];
-          for(int windowJ=-theDim/2;windowJ<(theDim+1)/2;++windowJ){
-            for(int windowI=-theDim/2;windowI<(theDim+1)/2;++windowI){
-              if(disc_opt[0]&&(windowI*windowI+windowJ*windowJ>(theDim/2)*(theDim/2)))
-                continue;
-              int j=j_centre+windowJ;
-              //check if j is out of bounds
-              if(static_cast<int>(j)<0||static_cast<int>(j)>=inputReader.nrOfRow())
-                continue;
-              int i=i_centre+windowI;
-              //check if i is out of bounds
-              if(static_cast<int>(i)<0||static_cast<int>(i)>=inputReader.nrOfCol())
-                continue;
-              if(verbose_opt[0])
-                cout << setprecision(12) << "reading image value at x,y " << x << "," << y << " (" << i << "," << j << "), ";
-              inputReader.readData(inputValue,GDT_Int16,i,j,band_opt[0]);
-              inputValues.push_back(inputValue);
-              if(inputValues.back()!=*(inputValues.begin()))
-                isHomogeneous=false;
-              if(verbose_opt[0])
-                cout << "input value: " << inputValue << endl;
-              pixelFlagged=false;
-              for(int iflag=0;iflag<nodata_opt.size();++iflag){
-                if(inputValue==nodata_opt[iflag]){
-                  pixelFlagged=true;
-                  break;
-                }
-              }
-              maskFlagged=false;//(masknodata_opt[ivalue]>=0)?false:true;
-              if(mask_opt.size()){
-                maskReader.readData(maskValue,GDT_Int16,i,j,band_opt[0]);
-                for(int ivalue=0;ivalue<masknodata_opt.size();++ivalue){
-                  if(masknodata_opt[ivalue]>=0){//values set in masknodata_opt are invalid
-                    if(maskValue==masknodata_opt[ivalue]){
-                      maskFlagged=true;
-                      break;
-                    }
-                  }
-                  else{//only values set in masknodata_opt are valid
-                    if(maskValue!=-masknodata_opt[ivalue])
-                      maskFlagged=true;
-                    else{
-                      maskFlagged=false;
-                      break;
-                    }
-                  }
-                }
-              }
-              pixelFlagged=pixelFlagged||maskFlagged;
-              if(pixelFlagged)
-                windowHasFlag=true;
-              else
-                windowAllFlagged=false;//at least one good pixel in neighborhood
-            }
-          }
-          //at this point we know the values for the entire window
-          if(homogeneous_opt[0]){//only centre pixel
-            int j=j_centre;
-            int i=i_centre;
-            //flag if not all pixels are homogeneous or if at least one pixel flagged
+	for(int ilayer=0;ilayer<referenceReaderOgr.getDataSource()->GetLayerCount();++ilayer){
+	  OGRLayer  *readLayer;
+	  readLayer = referenceReaderOgr.getDataSource()->GetLayer(ilayer);
+	  readLayer->ResetReading();
+	  OGRLayer *writeLayer;
+	  if(output_opt.size()){
+	    if(verbose_opt[0])
+	      cout << "creating output vector file " << output_opt[0] << endl;
+	    // assert(output_opt[0].find(".shp")!=string::npos);
+	    try{
+	      ogrWriter.open(output_opt[iref],ogrformat_opt[0]);
+	    }
+	    catch(string error){
+	      cerr << error << endl;
+	      exit(1);
+	    }
+	    char     **papszOptions=NULL;
+	    string layername=type2string<int>(ilayer);//output_opt[0].substr(0,output_opt[0].find(".shp"));
+	    if(verbose_opt[0])
+	      cout << "creating layer: " << layername << endl;
+	    if(ogrWriter.createLayer(layername, referenceReaderOgr.getProjection(ilayer), referenceReaderOgr.getGeometryType(ilayer), papszOptions)==NULL)
+	      cout << "Error: create layer failed!" << endl;
+	    else if(verbose_opt[0])
+	      cout << "created layer" << endl;
+	    if(verbose_opt[0])
+	      cout << "copy fields from " << reference_opt[iref] << endl;
+	    ogrWriter.copyFields(referenceReaderOgr,ilayer);
+	    //create extra field for classified label
+	    short theDim=boundary_opt[0];
+	    for(int windowJ=-theDim/2;windowJ<(theDim+1)/2;++windowJ){
+	      for(int windowI=-theDim/2;windowI<(theDim+1)/2;++windowI){
+		if(disc_opt[0]&&(windowI*windowI+windowJ*windowJ>(theDim/2)*(theDim/2)))
+		  continue;
+		ostringstream fs;
+		if(theDim>1)
+		  fs << labelclass_opt[0] << "_" << windowJ << "_" << windowI;
+		else
+		  fs << labelclass_opt[0];
+		if(verbose_opt[0])
+		  cout << "creating field " << fs.str() << endl;
+		ogrWriter.createField(fs.str(),OFTInteger,ilayer);
+	      }
+	    }
+	    writeLayer=ogrWriter.getDataSource()->GetLayer(ilayer);
+	  }
+	  OGRFeature *readFeature;
+	  int isample=0;
+	  while( (readFeature = readLayer->GetNextFeature()) != NULL ){
+	    if(verbose_opt[0])
+	      cout << "sample " << ++isample << endl;
+	    //get x and y from readFeature
+	    double x,y;
+	    OGRGeometry *poGeometry;
+	    OGRPoint centroidPoint;
+	    OGRPoint *poPoint;
+	    poGeometry = readFeature->GetGeometryRef();
+	    // assert( poGeometry != NULL && wkbFlatten(poGeometry->getGeometryType()) == wkbPoint );
+	    if(poGeometry==NULL)
+	      continue;
+	    else if(wkbFlatten(poGeometry->getGeometryType()) == wkbMultiPolygon){
+	      OGRMultiPolygon readPolygon = *((OGRMultiPolygon *) poGeometry);
+	      readPolygon = *((OGRMultiPolygon *) poGeometry);
+	      readPolygon.Centroid(&centroidPoint);
+	      poPoint=¢roidPoint;
+	    }
+	    else if(wkbFlatten(poGeometry->getGeometryType()) == wkbPolygon){
+	      OGRPolygon readPolygon=*((OGRPolygon *) poGeometry);
+	      readPolygon.Centroid(&centroidPoint);
+	      poPoint=¢roidPoint;
+	    }
+	    else if(wkbFlatten(poGeometry->getGeometryType()) == wkbPoint )
+	      poPoint = (OGRPoint *) poGeometry;
+	    else{
+	      std::cerr << "Warning: skipping feature (not of type point or polygon)" << std::endl;
+	      continue;
+	    }
+	    x=poPoint->getX();
+	    y=poPoint->getY();
+	    short inputValue;
+	    vector<short> inputValues;
+	    bool isHomogeneous=true;
+	    short maskValue;
+	    short outputValue;
+	    //read referenceValue from feature
+	    short referenceValue;
+	    string referenceClassName;
+	    if(classValueMap.size()){
+	      referenceClassName=readFeature->GetFieldAsString(readFeature->GetFieldIndex(labelref_opt[0].c_str()));
+	      referenceValue=classValueMap[referenceClassName];
+	    }
+	    else
+	      referenceValue=readFeature->GetFieldAsInteger(readFeature->GetFieldIndex(labelref_opt[0].c_str()));
+	    if(verbose_opt[0])
+	      cout << "reference value: " << referenceValue << endl;
+	    bool pixelFlagged=false;
+	    bool maskFlagged=false;
+	    for(int iflag=0;iflag<nodata_opt.size();++iflag){
+	      if(referenceValue==nodata_opt[iflag])
+		pixelFlagged=true;
+	    }
+	    if(pixelFlagged)
+	      continue;
+	    double i_centre,j_centre;
+	    //input reader is georeferenced!
+	    inputReader.geo2image(x,y,i_centre,j_centre);
+	    //       else{
+	    //         i_centre=x;
+	    //         j_centre=y;
+	    //       }
+	    //nearest neighbour
+	    j_centre=static_cast<int>(j_centre);
+	    i_centre=static_cast<int>(i_centre);
+	    //check if j_centre is out of bounds
+	    if(static_cast<int>(j_centre)<0||static_cast<int>(j_centre)>=inputReader.nrOfRow())
+	      continue;
+	    //check if i_centre is out of bounds
+	    if(static_cast<int>(i_centre)<0||static_cast<int>(i_centre)>=inputReader.nrOfCol())
+	      continue;
+	    OGRFeature *writeFeature;
+	    if(output_opt.size()){
+	      writeFeature = OGRFeature::CreateFeature(writeLayer->GetLayerDefn());
+	      if(verbose_opt[0])
+		cout << "copying fields from " << reference_opt[0] << endl;
+	      //todo: segmentation fault when more than one layer...?
+	      //test
+	      // assert(readFeature);
+	      // cout << readFeature->GetFieldAsString(readFeature->GetFieldIndex(labelref_opt[0].c_str())) << endl;
+	      // assert(writeFeature);
+	      // cout << "number of fields in writeFeature: " << writeFeature->GetFieldCount() << endl;
+	      //end test	      
+	      if(writeFeature->SetFrom(readFeature)!= OGRERR_NONE)
+		cerr << "writing feature failed" << endl;
+	      if(verbose_opt[0])
+		cout << "feature written" << endl;
+	    }
+	    bool windowAllFlagged=true;
+	    bool windowHasFlag=false;
+	    short theDim=boundary_opt[0];
+	    for(int windowJ=-theDim/2;windowJ<(theDim+1)/2;++windowJ){
+	      for(int windowI=-theDim/2;windowI<(theDim+1)/2;++windowI){
+		if(disc_opt[0]&&(windowI*windowI+windowJ*windowJ>(theDim/2)*(theDim/2)))
+		  continue;
+		int j=j_centre+windowJ;
+		//check if j is out of bounds
+		if(static_cast<int>(j)<0||static_cast<int>(j)>=inputReader.nrOfRow())
+		  continue;
+		int i=i_centre+windowI;
+		//check if i is out of bounds
+		if(static_cast<int>(i)<0||static_cast<int>(i)>=inputReader.nrOfCol())
+		  continue;
+		if(verbose_opt[0])
+		  cout << setprecision(12) << "reading image value at x,y " << x << "," << y << " (" << i << "," << j << "), ";
+		inputReader.readData(inputValue,GDT_Int16,i,j,band_opt[0]);
+		inputValues.push_back(inputValue);
+		if(inputValues.back()!=*(inputValues.begin()))
+		  isHomogeneous=false;
+		if(verbose_opt[0])
+		  cout << "input value: " << inputValue << endl;
+		pixelFlagged=false;
+		for(int iflag=0;iflag<nodata_opt.size();++iflag){
+		  if(inputValue==nodata_opt[iflag]){
+		    pixelFlagged=true;
+		    break;
+		  }
+		}
+		maskFlagged=false;//(masknodata_opt[ivalue]>=0)?false:true;
+		if(mask_opt.size()){
+		  maskReader.readData(maskValue,GDT_Int16,i,j,band_opt[0]);
+		  for(int ivalue=0;ivalue<masknodata_opt.size();++ivalue){
+		    if(masknodata_opt[ivalue]>=0){//values set in masknodata_opt are invalid
+		      if(maskValue==masknodata_opt[ivalue]){
+			maskFlagged=true;
+			break;
+		      }
+		    }
+		    else{//only values set in masknodata_opt are valid
+		      if(maskValue!=-masknodata_opt[ivalue])
+			maskFlagged=true;
+		      else{
+			maskFlagged=false;
+			break;
+		      }
+		    }
+		  }
+		}
+		pixelFlagged=pixelFlagged||maskFlagged;
+		if(pixelFlagged)
+		  windowHasFlag=true;
+		else
+		  windowAllFlagged=false;//at least one good pixel in neighborhood
+	      }
+	    }
+	    //at this point we know the values for the entire window
+	    if(homogeneous_opt[0]){//only centre pixel
+	      int j=j_centre;
+	      int i=i_centre;
+	      //flag if not all pixels are homogeneous or if at least one pixel flagged
           
-            if(!windowHasFlag&&isHomogeneous){
-              if(output_opt.size())
-                writeFeature->SetField(labelclass_opt[0].c_str(),static_cast<int>(inputValue));
-              if(confusion_opt[0]){
-                ++ntotalValidation;
-                if(classValueMap.size()){
-                  assert(inputValue<nameVector.size());
-                  string className=nameVector[inputValue];
-                  cm.incrementResult(type2string<short>(classValueMap[referenceClassName]),type2string<short>(classValueMap[className]),1);
-                }
-                else{
-                  int rc=distance(referenceRange.begin(),find(referenceRange.begin(),referenceRange.end(),referenceValue));
-                  int ic=distance(inputRange.begin(),find(inputRange.begin(),inputRange.end(),inputValue));
-                  assert(rc<nclass);
-                  assert(ic<nclass);
-                  ++nvalidation[rc];
-                  ++resultClass[rc][ic];
-                  if(verbose_opt[0]>1)
-                    cout << "increment: " << rc << " " << referenceRange[rc] << " " << ic << " " << inputRange[ic] << endl;
-                  cm.incrementResult(cm.getClass(rc),cm.getClass(ic),1);
-                }
-              }
-              if(inputValue==referenceValue){//correct
-		outputValue=valueE_opt[0];
-		if(nodata_opt.size()){
-		  if(valueE_opt[0]==nodata_opt[0])
-		    outputValue=inputValue;
+	      if(!windowHasFlag&&isHomogeneous){
+		if(output_opt.size())
+		  writeFeature->SetField(labelclass_opt[0].c_str(),static_cast<int>(inputValue));
+		if(confusion_opt[0]){
+		  ++ntotalValidation;
+		  if(classValueMap.size()){
+		    assert(inputValue<nameVector.size());
+		    string className=nameVector[inputValue];
+		    cm.incrementResult(type2string<short>(classValueMap[referenceClassName]),type2string<short>(classValueMap[className]),1);
+		  }
+		  else{
+		    int rc=distance(referenceRange.begin(),find(referenceRange.begin(),referenceRange.end(),referenceValue));
+		    int ic=distance(inputRange.begin(),find(inputRange.begin(),inputRange.end(),inputValue));
+		    assert(rc<nclass);
+		    assert(ic<nclass);
+		    ++nvalidation[rc];
+		    ++resultClass[rc][ic];
+		    if(verbose_opt[0]>1)
+		      cout << "increment: " << rc << " " << referenceRange[rc] << " " << ic << " " << inputRange[ic] << endl;
+		    cm.incrementResult(cm.getClass(rc),cm.getClass(ic),1);
+		  }
 		}
-              }
-              else if(inputValue>referenceValue)//1=forest,2=non-forest
-                outputValue=valueO_opt[0];//omission error
-              else
-                outputValue=valueC_opt[0];//commission error
-            }
-          }
-          else{
-            for(int windowJ=-theDim/2;windowJ<(theDim+1)/2;++windowJ){
-              for(int windowI=-theDim/2;windowI<(theDim+1)/2;++windowI){
-                if(disc_opt[0]&&(windowI*windowI+windowJ*windowJ>(theDim/2)*(theDim/2)))
-                  continue;
-                int j=j_centre+windowJ;
-                //check if j is out of bounds
-                if(static_cast<int>(j)<0||static_cast<int>(j)>=inputReader.nrOfRow())
-                  continue;
-                int i=i_centre+windowI;
-                //check if i is out of bounds
-                if(static_cast<int>(i)<0||static_cast<int>(i)>=inputReader.nrOfCol())
-                  continue;
-                if(!windowAllFlagged){
-                  ostringstream fs;
-                  if(theDim>1)
-                    fs << labelclass_opt[0] << "_" << windowJ << "_" << windowI;
-                  else
-                    fs << labelclass_opt[0];
-                  if(output_opt.size())
-                    writeFeature->SetField(fs.str().c_str(),static_cast<int>(inputValue));
-                  if(!windowJ&&!windowI){//centre pixel
-                    if(confusion_opt[0]){
-                      ++ntotalValidation;
-                      if(classValueMap.size()){
-                        assert(inputValue<nameVector.size());
-                        string className=nameVector[inputValue];
-                        cm.incrementResult(type2string<short>(classValueMap[referenceClassName]),type2string<short>(classValueMap[className]),1);
-                      }
-                      else{
-                        int rc=distance(referenceRange.begin(),find(referenceRange.begin(),referenceRange.end(),referenceValue));
-                        int ic=distance(inputRange.begin(),find(inputRange.begin(),inputRange.end(),inputValue));
-                        if(rc>=nclass)
-                          continue;
-                        if(ic>=nclass)
-                          continue;
-                        // assert(rc<nclass);
-                        // assert(ic<nclass);
-                        ++nvalidation[rc];
-                        ++resultClass[rc][ic];
-                        if(verbose_opt[0]>1)
-                          cout << "increment: " << rc << " " << referenceRange[rc] << " " << ic << " " << inputRange[ic] << endl;
-                        cm.incrementResult(cm.getClass(rc),cm.getClass(ic),1);
-                      }
-                    }
-                    if(inputValue==referenceValue){//correct
-		      outputValue=valueE_opt[0];
-		      if(nodata_opt.size()){
-			if(valueE_opt[0]==nodata_opt[0])
-			  outputValue=inputValue;
+		if(inputValue==referenceValue){//correct
+		  outputValue=valueE_opt[0];
+		  if(nodata_opt.size()){
+		    if(valueE_opt[0]==nodata_opt[0])
+		      outputValue=inputValue;
+		  }
+		}
+		else if(inputValue>referenceValue)//1=forest,2=non-forest
+		  outputValue=valueO_opt[0];//omission error
+		else
+		  outputValue=valueC_opt[0];//commission error
+	      }
+	    }
+	    else{
+	      for(int windowJ=-theDim/2;windowJ<(theDim+1)/2;++windowJ){
+		for(int windowI=-theDim/2;windowI<(theDim+1)/2;++windowI){
+		  if(disc_opt[0]&&(windowI*windowI+windowJ*windowJ>(theDim/2)*(theDim/2)))
+		    continue;
+		  int j=j_centre+windowJ;
+		  //check if j is out of bounds
+		  if(static_cast<int>(j)<0||static_cast<int>(j)>=inputReader.nrOfRow())
+		    continue;
+		  int i=i_centre+windowI;
+		  //check if i is out of bounds
+		  if(static_cast<int>(i)<0||static_cast<int>(i)>=inputReader.nrOfCol())
+		    continue;
+		  if(!windowAllFlagged){
+		    ostringstream fs;
+		    if(theDim>1)
+		      fs << labelclass_opt[0] << "_" << windowJ << "_" << windowI;
+		    else
+		      fs << labelclass_opt[0];
+		    if(output_opt.size())
+		      writeFeature->SetField(fs.str().c_str(),static_cast<int>(inputValue));
+		    if(!windowJ&&!windowI){//centre pixel
+		      if(confusion_opt[0]){
+			++ntotalValidation;
+			if(classValueMap.size()){
+			  assert(inputValue<nameVector.size());
+			  string className=nameVector[inputValue];
+			  cm.incrementResult(type2string<short>(classValueMap[referenceClassName]),type2string<short>(classValueMap[className]),1);
+			}
+			else{
+			  int rc=distance(referenceRange.begin(),find(referenceRange.begin(),referenceRange.end(),referenceValue));
+			  int ic=distance(inputRange.begin(),find(inputRange.begin(),inputRange.end(),inputValue));
+			  if(rc>=nclass)
+			    continue;
+			  if(ic>=nclass)
+			    continue;
+			  // assert(rc<nclass);
+			  // assert(ic<nclass);
+			  ++nvalidation[rc];
+			  ++resultClass[rc][ic];
+			  if(verbose_opt[0]>1)
+			    cout << "increment: " << rc << " " << referenceRange[rc] << " " << ic << " " << inputRange[ic] << endl;
+			  cm.incrementResult(cm.getClass(rc),cm.getClass(ic),1);
+			}
 		      }
-                    }
-                    else if(inputValue>referenceValue)//1=forest,2=non-forest
-                      outputValue=valueO_opt[0];//omission error
-                    else
-                      outputValue=valueC_opt[0];//commission error
-                  }
-                }
-              }
-            }
-          }
-          if(output_opt.size()){
-            if(!windowAllFlagged){
-              if(verbose_opt[0])
-                cout << "creating feature" << endl;
-              if(writeLayer->CreateFeature( writeFeature ) != OGRERR_NONE ){
-                string errorString="Failed to create feature in shapefile";
-                throw(errorString);
-              }
-            }
-            OGRFeature::DestroyFeature( writeFeature );
-          }
-        }
+		      if(inputValue==referenceValue){//correct
+			outputValue=valueE_opt[0];
+			if(nodata_opt.size()){
+			  if(valueE_opt[0]==nodata_opt[0])
+			    outputValue=inputValue;
+			}
+		      }
+		      else if(inputValue>referenceValue)//1=forest,2=non-forest
+			outputValue=valueO_opt[0];//omission error
+		      else
+			outputValue=valueC_opt[0];//commission error
+		    }
+		  }
+		}
+	      }
+	    }
+	    if(output_opt.size()){
+	      if(!windowAllFlagged){
+		if(verbose_opt[0])
+		  cout << "creating feature" << endl;
+		if(writeLayer->CreateFeature( writeFeature ) != OGRERR_NONE ){
+		  string errorString="Failed to create feature in shapefile";
+		  throw(errorString);
+		}
+	      }
+	      OGRFeature::DestroyFeature( writeFeature );
+	    }
+	  }
+	}
         if(output_opt.size())
           ogrWriter.close();
-        referenceReader.close();
+        referenceReaderOgr.close();
         inputReader.close();
         if(mask_opt.size())
           maskReader.close();
       }
     }
-  }//reference is shape file
+  }//reference is OGR
   else{
-    ImgWriterGdal imgWriter;
+    ImgWriterGdal gdalWriter;
     try{
       inputReader.open(input_opt[0]);//,imagicX_opt[0],imagicY_opt[0]);
       if(mask_opt.size())
@@ -528,18 +561,18 @@ int main(int argc, char *argv[])
           theInterleave+=inputReader.getInterleave();
           option_opt.push_back(theInterleave);
         }
-        imgWriter.open(output_opt[0],inputReader.nrOfCol(),inputReader.nrOfRow(),1,inputReader.getDataType(),inputReader.getImageType(),option_opt);
+        gdalWriter.open(output_opt[0],inputReader.nrOfCol(),inputReader.nrOfRow(),1,inputReader.getDataType(),inputReader.getImageType(),option_opt);
 	if(nodata_opt.size())
-	  imgWriter.GDALSetNoDataValue(nodata_opt[0]);
+	  gdalWriter.GDALSetNoDataValue(nodata_opt[0]);
         if(inputReader.isGeoRef()){
-          imgWriter.copyGeoTransform(inputReader);
+          gdalWriter.copyGeoTransform(inputReader);
         }
         if(colorTable_opt.size())
-          imgWriter.setColorTable(colorTable_opt[0]);
+          gdalWriter.setColorTable(colorTable_opt[0]);
         else if(inputReader.getColorTable()!=NULL){
           if(verbose_opt[0])
             cout << "set colortable from input image" << endl;
-          imgWriter.setColorTable(inputReader.getColorTable());
+          gdalWriter.setColorTable(inputReader.getColorTable());
         }
       }
       else if(verbose_opt[0])
@@ -560,22 +593,22 @@ int main(int argc, char *argv[])
     int irow=0;
     int icol=0;
     double oldreferencerow=-1;
-    ImgReaderGdal referenceReader;
+    ImgReaderGdal referenceReaderGdal;
     try{
-      referenceReader.open(reference_opt[0]);//,rmagicX_opt[0],rmagicY_opt[0]);
+      referenceReaderGdal.open(reference_opt[0]);//,rmagicX_opt[0],rmagicY_opt[0]);
     }
     catch(string error){
       cerr << error << endl;
       exit(1);
     }
     if(inputReader.isGeoRef()){
-      assert(referenceReader.isGeoRef());
-      if(inputReader.getProjection()!=referenceReader.getProjection())
+      assert(referenceReaderGdal.isGeoRef());
+      if(inputReader.getProjection()!=referenceReaderGdal.getProjection())
         cout << "projection of input image and reference image are different!" << endl;
     }
-    vector<short> lineReference(referenceReader.nrOfCol());
+    vector<short> lineReference(referenceReaderGdal.nrOfCol());
     if(confusion_opt[0]){
-      referenceReader.getRange(referenceRange,band_opt[0]);
+      referenceReaderGdal.getRange(referenceRange,band_opt[0]);
       for(int iflag=0;iflag<nodata_opt.size();++iflag){
         vector<short>::iterator fit;
         fit=find(referenceRange.begin(),referenceRange.end(),nodata_opt[iflag]);
@@ -609,22 +642,22 @@ int main(int argc, char *argv[])
       for(icol=0;icol<inputReader.nrOfCol();++icol){
         //find col in reference
         inputReader.image2geo(icol,irow,x,y);
-        referenceReader.geo2image(x,y,ireference,jreference);
-        if(ireference<0||ireference>=referenceReader.nrOfCol()){
+        referenceReaderGdal.geo2image(x,y,ireference,jreference);
+        if(ireference<0||ireference>=referenceReaderGdal.nrOfCol()){
           cerr << ireference << " out of reference range!" << endl;
           cerr << x << " " << y << " " << icol << " " << irow << endl;
           cerr << x << " " << y << " " << ireference << " " << jreference << endl;
           exit(1);
         }
         if(jreference!=oldreferencerow){
-          if(jreference<0||jreference>=referenceReader.nrOfRow()){
+          if(jreference<0||jreference>=referenceReaderGdal.nrOfRow()){
             cerr << jreference << " out of reference range!" << endl;
             cerr << x << " " << y << " " << icol << " " << irow << endl;
             cerr << x << " " << y << " " << ireference << " " << jreference << endl;
             exit(1);
           }
           else{
-            referenceReader.readData(lineReference,GDT_Int16,static_cast<int>(jreference),band_opt[0]);
+            referenceReaderGdal.readData(lineReference,GDT_Int16,static_cast<int>(jreference),band_opt[0]);
             oldreferencerow=jreference;
           }
         }
@@ -703,11 +736,11 @@ int main(int argc, char *argv[])
       }
       if(output_opt.size()){
         try{
-          imgWriter.writeData(lineOutput,GDT_Int16,irow);
+          gdalWriter.writeData(lineOutput,GDT_Int16,irow);
         }
         catch(string errorstring){
           cerr << "lineOutput.size(): " << lineOutput.size() << endl;
-          cerr << "imgWriter.nrOfCol(): " << imgWriter.nrOfCol() << endl;
+          cerr << "gdalWriter.nrOfCol(): " << gdalWriter.nrOfCol() << endl;
           cerr << errorstring << endl;
           exit(1);
         }
@@ -722,14 +755,14 @@ int main(int argc, char *argv[])
         pfnProgress(progress,pszMessage,pProgressArg);
     }
     if(output_opt.size())
-      imgWriter.close();
+      gdalWriter.close();
     else if(!confusion_opt[0]){
       if(isDifferent)
         cout << input_opt[0] << " and " << reference_opt[0] << " are different" << endl;
       else
         cout << input_opt[0] << " and " << reference_opt[0] << " are identical" << endl;
     }
-    referenceReader.close();
+    referenceReaderGdal.close();
     inputReader.close();
     if(mask_opt.size())
       maskReader.close();
diff --git a/src/apps/pkeditogr.cc b/src/apps/pkeditogr.cc
index 44158ac..8de2f3a 100644
--- a/src/apps/pkeditogr.cc
+++ b/src/apps/pkeditogr.cc
@@ -33,9 +33,10 @@ int main(int argc, char *argv[])
   Optionpk<string> output_opt("o", "output", "Output vector file");
   Optionpk<string> ogrformat_opt("f", "f", "Output OGR file format","ESRI Shapefile");
   Optionpk<string> selectField_opt("select", "select", "select field (combined with like opt)");
+  //selectField can also be done via ogr2ogr using the -select option
   Optionpk<string> like_opt("like", "like", "substring(s) to be found in select field. If multiple substrings are provided, feature will be selected if one of them is found (stringent option must be false)");
   Optionpk<bool> stringent_opt("st", "stringent", "string in like option must exactly match to select feature)",false);
-  Optionpk<string> field_opt("f", "field", "output field names (number must exactly match input fields)");
+  // Optionpk<string> field_opt("as_field", "_field", "output field names (number must exactly match input fields)");
   //renaming fields can also be done via ogr2ogr using the -sql option:
   //ogr2ogr outdataset indataset -sql "SELECT src_field1 AS dst_field1, src_field2 AS dst_field2 FROM sourcelayer"
   Optionpk<long int> setfeature_opt("sf", "sf", "id of feature(s) to set (start from 0)");
@@ -54,7 +55,7 @@ int main(int argc, char *argv[])
     selectField_opt.retrieveOption(argc,argv);
     like_opt.retrieveOption(argc,argv);
     stringent_opt.retrieveOption(argc,argv);
-    field_opt.retrieveOption(argc,argv);
+    // field_opt.retrieveOption(argc,argv);
     addname_opt.retrieveOption(argc,argv);
     addtype_opt.retrieveOption(argc,argv);
     addvalue_opt.retrieveOption(argc,argv);
@@ -102,8 +103,8 @@ int main(int argc, char *argv[])
   if(verbose_opt[0])
     cout << "reset reading" << endl;
   ogrReader.getLayer()->ResetReading();
-  if(field_opt.size())
-    assert(field_opt.size()==ogrReader.getFieldCount());
+  // if(field_opt.size())
+  //   assert(field_opt.size()==ogrReader.getFieldCount());
   unsigned long int ifeature=0;
   if(verbose_opt[0])
     cout << "going through features" << endl << flush;
@@ -116,16 +117,16 @@ int main(int argc, char *argv[])
   writeFields=readFields;
   try{
     for(int ifield=0;ifield<readFields.size();++ifield){
-      if(field_opt.size()>ifield)
-        writeFields[ifield]->SetName(field_opt[ifield].c_str());
+      // if(field_opt.size()>ifield)
+      //   writeFields[ifield]->SetName(field_opt[ifield].c_str());
       if(verbose_opt[0])
         std::cout << readFields[ifield]->GetNameRef() << " -> " << writeFields[ifield]->GetNameRef() << std::endl;
       if(writeLayer->CreateField(writeFields[ifield]) != OGRERR_NONE ){
         ostringstream es;
-        if(field_opt.size()>ifield)
-          es << "Creating field " << field_opt[ifield] << " failed";
-        else
-          es << "Creating field " << readFields[ifield] << " failed";
+        // if(field_opt.size()>ifield)
+        //   es << "Creating field " << field_opt[ifield] << " failed";
+        // else
+	es << "Creating field " << readFields[ifield] << " failed";
         string errorString=es.str();
         throw(errorString);
       }
diff --git a/src/apps/pkreclass.cc b/src/apps/pkreclass.cc
index 2247a2a..21fefd8 100644
--- a/src/apps/pkreclass.cc
+++ b/src/apps/pkreclass.cc
@@ -115,10 +115,19 @@ int main(int argc, char *argv[])
     for(mit=codemapString.begin();mit!=codemapString.end();++mit)
       cout << (*mit).first << " " << (*mit).second << endl;
   }
-  if(input_opt[0].find(".shp")!=string::npos){//shape file
+  bool refIsRaster=false;
+  ImgReaderOgr ogrReader;
+  try{
+    ogrReader.open(input_opt[0]);
+  }
+  catch(string errorString){
+    refIsRaster=true;
+  }
+  // if(input_opt[0].find(".shp")!=string::npos){//shape file
+  if(!refIsRaster){
     if(verbose_opt[0])
       cout << "opening " << input_opt[0] << " for reading " << endl;
-    ImgReaderOgr ogrReader(input_opt[0]);
+    // ImgReaderOgr ogrReader(input_opt[0]);
     if(verbose_opt[0])
       cout << "opening " << output_opt[0] << " for writing " << endl;
     ImgWriterOgr ogrWriter(output_opt[0],ogrReader);
diff --git a/src/imageclasses/ImgWriterOgr.cc b/src/imageclasses/ImgWriterOgr.cc
index 15b769c..357472c 100644
--- a/src/imageclasses/ImgWriterOgr.cc
+++ b/src/imageclasses/ImgWriterOgr.cc
@@ -250,7 +250,7 @@ void ImgWriterOgr::copyFields(const ImgReaderOgr& imgReaderOgr, int theLayer){
   //get fields from imgReaderOgr
   std::vector<OGRFieldDefn*> fields;
   
-  imgReaderOgr.getFields(fields);
+  imgReaderOgr.getFields(fields,theLayer);
 //   OGRFeatureDefn *poFDefn = poLayer->GetLayerDefn();
   for(int iField=0;iField<fields.size();++iField){
     if(m_datasource->GetLayer(theLayer)->CreateField(fields[iField]) != OGRERR_NONE ){

-- 
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