#include <sstream>
#include <iostream>

#include "core/image/FilterT.h"
#include "core/basics/numerictools.h"
#include "core/basics/StringTools.h"
#include "core/basics/Timer.h"

#include <vislearning/classifier/fpclassifier/gphik/FPCGPHIK.h>
#include "vislearning/baselib/ICETools.h"
#include "vislearning/baselib/Globals.h"
#include "vislearning/features/fpfeatures/SparseVectorFeature.h"

#include "segmentation/GenericRegionSegmentationMethodSelection.h"

#include "SemSegNovelty.h"

using namespace std;
using namespace NICE;
using namespace OBJREC;

SemSegNovelty::SemSegNovelty ( )
    : SemanticSegmentation ( )
{
  this->forbidden_classesTrain.clear();
  this->forbidden_classesActiveLearning.clear();
  this->classesInUse.clear();
  
  this->globalMaxUncert = -numeric_limits<double>::max();  
  
  //we don't have queried any region so far
  this->queriedRegions.clear(); 
  
  this->featExtract = new LocalFeatureColorWeijer ();
  
  // those two guys need to be NULL, since only one of them will be active later on
  this->classifier = NULL;
  this->vclassifier = NULL;
  
  // this one here as well
  this->regionSeg = NULL;
}

SemSegNovelty::SemSegNovelty ( const Config * _conf,
                               const MultiDataset *md )
{
  ///////////
  // same code as in empty constructor - duplication can be avoided with C++11 allowing for constructor delegation
  ///////////  
  this->forbidden_classesTrain.clear();
  this->forbidden_classesActiveLearning.clear();
  this->classesInUse.clear();
  
  this->globalMaxUncert = -numeric_limits<double>::max();  
  
  //we don't have queried any region so far
  this->queriedRegions.clear(); 
  
  this->featExtract = new LocalFeatureColorWeijer ();
  
  // those two guys need to be NULL, since only one of them will be active later on
  this->classifier = NULL;
  this->vclassifier = NULL;
  
  // this one here as well
  this->regionSeg = NULL;  
  
  ///////////
  // here comes the new code part different from the empty constructor
  ///////////    
  this->setClassNames ( & ( md->getClassNames ( "train" ) ) );
  
  this->initFromConfig( _conf ); 
}

SemSegNovelty::~SemSegNovelty()
{
  if(newTrainExamples.size() > 0)
  {
    // show most uncertain region
    if (b_visualizeALimages)
      showImage(maskedImg);
    
    //incorporate new information into the classifier
    if (classifier != NULL)
    {
      //NOTE dangerous!
      classifier->addMultipleExamples(newTrainExamples);
    }
    
    //store the classifier, such that we can read it again in the next round (if we like that)
    classifier->save ( cache + "/classifier.data" );
  }
  
  // clean-up
  
  ///////////////////////////////
  //     FEATURE EXTRACTION    //
  ///////////////////////////////   
  if ( featExtract != NULL )
    delete featExtract;

  ///////////////////////////////
  //     CLASSIFICATION STUFF  //
  ///////////////////////////////   
  if ( classifier != NULL )
    delete classifier;
  if ( vclassifier != NULL )
    delete vclassifier;
  
  ///////////////////////////////
  //     SEGMENTATION STUFF    //
  ///////////////////////////////   
  if ( this->regionSeg != NULL )
    delete this->regionSeg;
}

void SemSegNovelty::initFromConfig(const Config* conf, const string _confSection) 
{
  //first of all, call method of parent object
  SemanticSegmentation::initFromConfig( conf );
  
  featExtract->initFromConfig ( conf );

  //save and read segmentation results from files
  this->reuseSegmentation = conf->gB ( "FPCPixel", "reuseSegmentation", true );
  //save the classifier to a file
  this->save_classifier = conf->gB ( "FPCPixel", "save_classifier", true ); 
  //read the classifier from a file
  this->read_classifier = conf->gB ( "FPCPixel", "read_classifier", false ); 

  //write uncertainty results in the same folder as done for the segmentation results
  resultdir = conf->gS("debug", "resultdir", "result");
  cache = conf->gS ( "cache", "root", "" );
  cache = conf->getAbsoluteFilenameRelativeToThisConfig(cache);

  this->findMaximumUncert = conf->gB(_confSection, "findMaximumUncert", true);
  this->whs = conf->gI ( _confSection, "window_size", 10 );
  //distance to next descriptor during training
  this->trainWSize = conf->gI ( _confSection, "train_window_size", 10 );
  //distance to next descriptor during testing
  this->testWSize = conf->gI (_confSection, "test_window_size", 10);
  // select your segmentation method here
  this->s_rsMethode = conf->gS ( _confSection, "segmentation", "none" );
 
  if( this->s_rsMethode == "none" )
  {
    regionSeg = NULL;
  }
  else
  {
    RegionSegmentationMethod *tmpRegionSeg = GenericRegionSegmentationMethodSelection::selectRegionSegmentationMethod( conf, this->s_rsMethode );    
    if ( reuseSegmentation )
      regionSeg = new RSCache ( conf, tmpRegionSeg );
    else
      regionSeg = tmpRegionSeg;
  }
  
  //define which measure for "novelty" we want to use
  noveltyMethodString = conf->gS( _confSection,  "noveltyMethod", "gp-variance");
  if (noveltyMethodString.compare("gp-variance") == 0)  // novel = large variance
  {
    this->noveltyMethod = GPVARIANCE;
    this->mostNoveltyWithMaxScores = true;
  }
  else if (noveltyMethodString.compare("gp-uncertainty") == 0) //novel = large uncertainty (mean / var)
  {
    this->noveltyMethod = GPUNCERTAINTY;
    this->mostNoveltyWithMaxScores = false;
    globalMaxUncert = numeric_limits<double>::max();
  } 
  else if (noveltyMethodString.compare("gp-mean") == 0) //novel = small mean
  {
    this->noveltyMethod = GPMINMEAN;
    this->mostNoveltyWithMaxScores = false;
    globalMaxUncert = numeric_limits<double>::max();
  }
  else if (noveltyMethodString.compare("gp-meanRatio") == 0)  //novel = small difference between mean of most plausible class and mean of snd
                                                              //        most plausible class (not useful in binary settings)
  {
    this->noveltyMethod = GPMEANRATIO;
    this->mostNoveltyWithMaxScores = false;
    globalMaxUncert = numeric_limits<double>::max();
  }
  else if (noveltyMethodString.compare("gp-weightAll") == 0) // novel = large weight in alpha vector after updating the model (can be predicted exactly)
  {
    this->noveltyMethod = GPWEIGHTALL;
    this->mostNoveltyWithMaxScores = true;
  }
  else if (noveltyMethodString.compare("gp-weightRatio") == 0) // novel = small difference between weights for alpha vectors 
                                                               //     with assumptions of GT label to be the most 
                                                               //     plausible against the second most plausible class   
  {
    this->noveltyMethod = GPWEIGHTRATIO;
    this->mostNoveltyWithMaxScores = false;
    globalMaxUncert = numeric_limits<double>::max();
  }
  else if (noveltyMethodString.compare("random") == 0) 
  {
     initRand(); 
     this->noveltyMethod = RANDOM;
  }
  else
  {
    this->noveltyMethod = GPVARIANCE;
    this->mostNoveltyWithMaxScores = true;
  }
  
  b_visualizeALimages = conf->gB(_confSection, "visualizeALimages", false);
  
  
  classifierString = conf->gS ( _confSection, "classifier", "GPHIKClassifier" );  
  classifier = NULL;
  vclassifier = NULL;
  if ( classifierString.compare("GPHIKClassifier") == 0)
  { 
    //just to make sure, that we do NOT perform an optimization after every iteration step
    //this would just take a lot of time, which is not desired so far
    //TODO edit this!
    //this->conf->sB( "GPHIKClassifier", "performOptimizationAfterIncrement", false );    
    classifier = new FPCGPHIK ( conf, "GPHIKClassifier" );
  }
  else
    vclassifier = GenericClassifierSelection::selectVecClassifier ( conf, classifierString ); 
  
  
  
  
  //check the same thing for the training classes - this is very specific to our setup 
  std::string forbidden_classesTrain_s = conf->gS ( "analysis", "donttrainTrain", "" );
  if ( forbidden_classesTrain_s == "" )
  {
    forbidden_classesTrain_s = conf->gS ( "analysis", "forbidden_classesTrain", "" );
  }
  this->classNames->getSelection ( forbidden_classesTrain_s, forbidden_classesTrain );  
}



void SemSegNovelty::visualizeRegion(const NICE::ColorImage &img, const NICE::Matrix &regions, int region, NICE::ColorImage &outimage)
{
  std::vector<uchar> color;
  color.push_back(255);
  color.push_back(0);
  color.push_back(0);
    
  int width = img.width();
  int height = img.height();
  
  outimage.resize(width,height);
  
  for(int y = 0; y < height; y++)
  {
    for(int x = 0; x < width; x++)
    {
      if(regions(x,y) == region)
      {
        for(int c = 0; c < 3; c++)
        {
          outimage(x,y,c) = color[c];
        }
      }
      else
      {
        for(int c = 0; c < 3; c++)
        {
          outimage(x,y,c) = img(x,y,c);
        }
      }
    }
  }
}

void SemSegNovelty::train ( const MultiDataset *md )
{
  if ( this->read_classifier )
  {
    try
    {
      if ( this->classifier != NULL )
      {
        string classifierdst = "/classifier.data";        
        fprintf ( stderr, "SemSegNovelty:: Reading classifier data from %s\n", ( cache + classifierdst ).c_str() );        
        classifier->read ( cache + classifierdst );
      }
      else
      {
        string classifierdst = "/veccl.data";        
        fprintf ( stderr, "SemSegNovelty:: Reading classifier data from %s\n", ( cache + classifierdst ).c_str() );          
        vclassifier->read ( cache + classifierdst );      
      }
      

      fprintf ( stderr, "SemSegNovelty:: successfully read\n" );
    }
    catch ( char *str )
    {
      cerr << "error reading data: " << str << endl;
    }
  }
  else
  {
    const LabeledSet train = * ( *md ) ["train"];
    const LabeledSet *trainp = &train;

    ////////////////////////
    // feature extraction //
    ////////////////////////
  
    ProgressBar pb ( "Local Feature Extraction" );
    pb.show();

    int imgnb = 0;

    Examples examples;
    examples.filename = "training";

    int featdim = -1;

    classesInUse.clear();  
    
    LOOP_ALL_S ( *trainp )
    {
      //EACH_S(classno, currentFile);
      EACH_INFO ( classno, info );

      std::string currentFile = info.img();

      CachedExample *ce = new CachedExample ( currentFile );

      const LocalizationResult *locResult = info.localization();
      if ( locResult->size() <= 0 )
      {
	fprintf ( stderr, "WARNING: NO ground truth polygons found for %s !\n",
		  currentFile.c_str() );
	continue;
      }

      int xsize, ysize;
      ce->getImageSize ( xsize, ysize );

      Image labels ( xsize, ysize );
      labels.set ( 0 );
      locResult->calcLabeledImage ( labels, ( *classNames ).getBackgroundClass() );

      NICE::ColorImage img;
      try {
	img = ColorImage ( currentFile );
      } catch ( Exception ) {
	cerr << "SemSegNovelty: error opening image file <" << currentFile << ">" << endl;
	continue;
      }

      Globals::setCurrentImgFN ( currentFile );

      MultiChannelImageT<double> feats;

      // extract features
      featExtract->getFeats ( img, feats );
      featdim = feats.channels();
      feats.addChannel(featdim);

      for (int c = 0; c < featdim; c++)
      {
	ImageT<double> tmp = feats[c];
	ImageT<double> tmp2 = feats[c+featdim];

	NICE::FilterT<double, double, double>::gradientStrength (tmp, tmp2);
      }
      featdim += featdim;

      // compute integral images
      for ( int c = 0; c < featdim; c++ )
      {
	feats.calcIntegral ( c );
      }

      for ( int y = 0; y < ysize; y += trainWSize)
      {
    for ( int x = 0; x < xsize; x += trainWSize )
	{

	  int classnoTmp = labels.getPixel ( x, y );
	  
	  if ( forbidden_classesTrain.find ( classnoTmp ) != forbidden_classesTrain.end() )
	  {
	    continue;
	  }
	  
	  if (classesInUse.find(classnoTmp) == classesInUse.end())
	  {
	    classesInUse.insert(classnoTmp);
	  }
	  
	  Example example;
	  example.vec = NULL;
	  example.svec = new SparseVector ( featdim );
	  for ( int f = 0; f < featdim; f++ )
	  {
	    double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
	    if ( val > 1e-10 )
	      ( *example.svec ) [f] = val;
	  }

	  example.svec->normalize();

	  example.position = imgnb;
	  examples.push_back ( pair<int, Example> ( classnoTmp, example ) );

	}
      }
  
      
      

      delete ce;
      imgnb++;
      pb.update ( trainp->count() );
    }
    
      
    numberOfClasses = classesInUse.size();
    std::cerr << "numberOfClasses: " << numberOfClasses << std::endl;  
    std::cerr << "classes in use: " << std::endl;
    for (std::set<int>::const_iterator it = classesInUse.begin(); it != classesInUse.end(); it++)
    {
      std::cerr << *it << " ";
    }    
    std::cerr << std::endl;

    pb.hide();


    //////////////////////
    // train classifier //
    //////////////////////
    FeaturePool fp;

    Feature *f = new SparseVectorFeature ( featdim );

    f->explode ( fp );
    delete f;

    if ( classifier != NULL )
    {
      std::cerr << "train FP-classifier with " << examples.size() << " examples" << std::endl;
      classifier->train ( fp, examples );
      std::cerr << "training finished" << std::endl;
    }
    else
    {
      LabeledSetVector lvec;
      convertExamplesToLSet ( examples, lvec );
      vclassifier->teach ( lvec );
  //     if ( usegmm )
  //       convertLSetToSparseExamples ( examples, lvec );
  //     else
      std::cerr << "classifierString: " << classifierString << std::endl;
      if (this->classifierString.compare("nn") == 0)
      {
	convertLSetToExamples ( examples, lvec, true /* only remove pointers to the data in the LSet-struct*/);
      }
      else
      {
	convertLSetToExamples ( examples, lvec, false /* remove all training examples of the LSet-struct */);
      }
      vclassifier->finishTeaching();
    }  

    fp.destroy();

    if ( save_classifier )
    {
      if ( classifier != NULL )
	classifier->save ( cache + "/classifier.data" );
      else
	vclassifier->save ( cache + "/veccl.data" );    
    }

    ////////////
    //clean up//
    ////////////
    for ( int i = 0; i < ( int ) examples.size(); i++ )
    {
      examples[i].second.clean();
    }
    examples.clear();

    cerr << "SemSeg training finished" << endl;
  }
}


void SemSegNovelty::semanticseg (
        CachedExample *ce,
        NICE::ImageT<int> & segresult,
        NICE::MultiChannelImageT<double> & probabilities )
{  
  Timer timer;
  timer.start();
  
  //segResult contains the GT labels when this method is called
  // we simply store them in labels, to have an easy access to the GT information lateron
  NICE::ImageT<int> labels = segresult;
  //just to be sure that we do not have a GT-biased result :)
  segresult.set(0);

  int featdim = -1;

  std::string currentFile = Globals::getCurrentImgFN();


  int xsize, ysize;
  ce->getImageSize ( xsize, ysize );

  probabilities.reInit( xsize, ysize, this->classNames->getMaxClassno() + 1);
  probabilities.setAll ( -std::numeric_limits<double>::max() );
   
  NICE::ColorImage img;
  try {
    img = ColorImage ( currentFile );
  } catch ( Exception ) {
    cerr << "SemSegNovelty: error opening image file <" << currentFile << ">" << endl;
    return;
  }

 // MultiChannelImageT<double> m_CurrentImageFeatures;

  // extract features
  m_CurrentImageFeatures.freeData();
  featExtract->getFeats ( img, m_CurrentImageFeatures );
  featdim = m_CurrentImageFeatures.channels();
  m_CurrentImageFeatures.addChannel(featdim);

  for (int c = 0; c < featdim; c++)
  {
    ImageT<double> tmp = m_CurrentImageFeatures[c];
    ImageT<double> tmp2 = m_CurrentImageFeatures[c+featdim];

    NICE::FilterT<double, double, double>::gradientStrength (tmp, tmp2);
  }
  featdim += featdim;

  // compute integral images
  for ( int c = 0; c < featdim; c++ )
  {
    m_CurrentImageFeatures.calcIntegral ( c );
  }
  
  timer.stop();
  std::cout << "AL time for preparation: " << timer.getLastAbsolute() << std::endl;
    
  timer.start();
  //classification results currently only needed to be computed separately if we use the vclassifier, i.e., the nearest neighbor used 
  // for the "novel feature learning" approach
  //in all other settings, such as active sem seg in general, we do this within the novelty-computation-methods
  if ( classifier == NULL )
  {
    this->computeClassificationResults( m_CurrentImageFeatures, segresult, probabilities, xsize, ysize, featdim);
  }
//   timer.stop();
//   
//   std::cerr << "classification results computed" << std::endl;
  
  FloatImage noveltyImage ( xsize, ysize );
  noveltyImage.set ( 0.0 );  
  
  switch (noveltyMethod)
  {
    case GPVARIANCE:
    {
         this->computeNoveltyByVariance( noveltyImage, m_CurrentImageFeatures, segresult, probabilities, xsize, ysize,  featdim );
         break;
    }
    case GPUNCERTAINTY:
    {
         this->computeNoveltyByGPUncertainty( noveltyImage, m_CurrentImageFeatures, segresult, probabilities, xsize, ysize,  featdim );
         break;         
    }
    case GPMINMEAN:
    {
         std::cerr << "compute novelty using the minimum mean" << std::endl;
         this->computeNoveltyByGPMean( noveltyImage, m_CurrentImageFeatures, segresult, probabilities, xsize, ysize,  featdim );
         break;         
    }
    case GPMEANRATIO:
    {
         this->computeNoveltyByGPMeanRatio( noveltyImage, m_CurrentImageFeatures, segresult, probabilities, xsize, ysize,  featdim );
         break;         
    }
    case GPWEIGHTALL:
    {
         this->computeNoveltyByGPWeightAll( noveltyImage, m_CurrentImageFeatures, segresult, probabilities, xsize, ysize,  featdim );
         break;         
    }
    case GPWEIGHTRATIO:
    {
         this->computeNoveltyByGPWeightRatio( noveltyImage, m_CurrentImageFeatures, segresult, probabilities, xsize, ysize,  featdim );
         break;         
    }    
    case RANDOM:
    {
         this->computeNoveltyByRandom( noveltyImage, m_CurrentImageFeatures, segresult, probabilities, xsize, ysize,  featdim );
         break;               
    }
    default:
    {
         //do nothing, keep the image constant to 0.0
         break;
    }
         
  }
  
  timer.stop();
  std::cout << "AL time for novelty score computation: " << timer.getLastAbsolute() << std::endl;
  
  if (b_visualizeALimages)
  {
      ColorImage imgrgbTmp (xsize, ysize);
      ICETools::convertToRGB ( noveltyImage, imgrgbTmp );
      showImage(imgrgbTmp, "Novelty Image without Region Segmentation");  
  }
  
    
  timer.start();
  
  //Regionen ermitteln
  if(regionSeg != NULL)
  {
    NICE::Matrix mask;
    int amountRegions = regionSeg->segRegions ( img, mask );
    
    //compute probs per region
    std::vector<std::vector<double> > regionProb(amountRegions, std::vector<double>(probabilities.channels(), -std::numeric_limits<double>::max() ));
    std::vector<double> regionNoveltyMeasure (amountRegions, 0.0);

    std::vector<int> regionCounter(amountRegions, 0);
    std::vector<int> regionCounterNovelty(amountRegions, 0);
    for ( int y = 0; y < ysize; y += testWSize) //y++)
    {
      for (int x = 0; x < xsize; x += testWSize) //x++)
      {
        int r = mask(x,y);
        regionCounter[r]++;

        for(int j = 0; j < probabilities.channels(); j++)
        {
            if( regionProb[r][j] == -std::numeric_limits<double>::max() )
                regionProb[r][j] = 0.0f;
            regionProb[r][j] += probabilities ( x, y, j );
        }
        
        if ( forbidden_classesActiveLearning.find( labels(x,y) ) == forbidden_classesActiveLearning.end() )
        {
          //count the amount of "novelty" for the corresponding region
          regionNoveltyMeasure[r] += noveltyImage(x,y);
          regionCounterNovelty[r]++;
        }
      }
    }
       
    //find best class per region
    std::vector<int> bestClassPerRegion(amountRegions,0);
    
    double maxNoveltyScore = -numeric_limits<double>::max();
    if (!mostNoveltyWithMaxScores)
    {
      maxNoveltyScore = numeric_limits<double>::max();
    }   
    
    int maxUncertRegion = -1;
    
    //loop over all regions and compute averaged novelty scores
    for(int r = 0; r < amountRegions; r++)
    {
      
      //check for the most plausible class per region
      double maxval = -numeric_limits<double>::max();
      
      //loop over all classes
      for(int c = 0; c < probabilities.channels(); c++)
      {
        regionProb[r][c] /= regionCounter[r];
        
        if(  (maxval < regionProb[r][c]) ) //&& (regionProb[r][c] != 0.0) ) 
        {        
              maxval = regionProb[r][c];
              bestClassPerRegion[r] = c;
        }
      }
       
      //if the region only contains unvalid information (e.g., background) skip it
      if (regionCounterNovelty[r] == 0)
      {
        continue;
      }
      
      //normalize summed novelty scores to region size
      regionNoveltyMeasure[r] /= regionCounterNovelty[r];
    
      //did we find a region that has a higher score as the most novel region known so far within this image?
      if(   (  mostNoveltyWithMaxScores && (maxNoveltyScore < regionNoveltyMeasure[r]) )    // if we look for large novelty scores, e.g., variance
        || ( !mostNoveltyWithMaxScores && (maxNoveltyScore > regionNoveltyMeasure[r]) ) )  // if we look for small novelty scores, e.g., min mean
      {
                   //did we already query a region of this image?                --   and it was this specific region
        if ( (queriedRegions.find( currentFile ) != queriedRegions.end() ) && ( queriedRegions[currentFile].find(r) != queriedRegions[currentFile].end() ) )
        {
          continue;
        }
        else //only accept the region as novel if we never queried it before
        {
          maxNoveltyScore = regionNoveltyMeasure[r];
          maxUncertRegion = r;        
        }

      }

    }
    
    // after finding the most novel region for the current image, check whether this region is also the most novel with respect
    // to all previously seen test images
    // if so, store the corresponding features, since we want to "actively" query them to incorporate useful information
    if(findMaximumUncert)
    {
      if(    (   mostNoveltyWithMaxScores && (maxNoveltyScore > globalMaxUncert) )
          || (  !mostNoveltyWithMaxScores && (maxNoveltyScore < globalMaxUncert) ) )
      {
        //current most novel region of the image has "higher" novelty score then previous most novel region of all test images worked on so far
        // -> save new important features of this region
        Examples examples;

        for ( int y = 0; y < ysize; y += testWSize )
        {
          for ( int x = 0; x < xsize; x += testWSize)
          {
              if(mask(x,y) == maxUncertRegion)
              {
                  int classnoTmp = labels(x,y);
                  if ( forbidden_classesActiveLearning.find(classnoTmp) != forbidden_classesActiveLearning.end() )
                      continue;
                  Example t_Example(NULL, x, y);
                  t_Example.vec = NULL;
                  t_Example.svec = new SparseVector ( featdim );
                  for ( int f = 0; f < featdim; f++ )
                  {
//                      double val =  ( *example.svec ) [f];
                       double val = m_CurrentImageFeatures.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
                    if ( val > 1e-10 )
                      ( *t_Example.svec ) [f] = val;
                  }
                  examples.push_back ( pair<int, Example> ( classnoTmp, t_Example ) );
              }
          }
        }
        
        if(examples.size() > 0)
        {
          std::cerr << "found " << examples.size() << " new examples in the queried region" << std::endl << std::endl;
          // sauber aufräumen
          for( int i=0; i< newTrainExamples.size(); i++)
          {
              delete newTrainExamples.at(i).second.svec;
              newTrainExamples.at(i).second.svec = NULL;
          }
          newTrainExamples.clear();
          newTrainExamples = examples;
          globalMaxUncert = maxNoveltyScore;
          //prepare for later visualization
//           if (b_visualizeALimages)
            visualizeRegion(img,mask,maxUncertRegion,maskedImg);
        }
        else
        {
          std::cerr << "the queried region has no valid information" << std::endl << std::endl;
        }
        
        //save filename and region index
        currentRegionToQuery.first = currentFile;
        currentRegionToQuery.second = maxUncertRegion;
      }
    }

    //write back best results per region
    //i.e., write normalized novelty scores for every region into the novelty image
    for ( int y = 0; y < ysize; y++)
    {
      for (int x = 0; x < xsize; x++)
      {
        int r = mask(x,y);
        for(int j = 0; j < probabilities.channels(); j++)
        {
          probabilities ( x, y, j ) = regionProb[r][j];
        }
        segresult(x,y) = bestClassPerRegion[r];
        // write novelty scores for every segment into the "final" image
        noveltyImage(x,y) = regionNoveltyMeasure[r];
      }
    }
  } // if regionSeg != null
  
  timer.stop();
  std::cout << "AL time for determination of novel regions: " << timer.getLastAbsolute() << std::endl;

  if (b_visualizeALimages)
  {
    timer.start();

    std::stringstream out;
    std::vector< std::string > list2;
    StringTools::split ( Globals::getCurrentImgFN (), '/', list2 );
    out << resultdir << "/" << list2.back();

    noveltyImage.writeRaw(out.str() + "_run_" +  NICE::intToString(this->iterationCountSuffix) + "_" + noveltyMethodString+".rawfloat");

    ColorImage imgrgb ( xsize, ysize );
    ICETools::convertToRGB ( noveltyImage, imgrgb );
    showImage(imgrgb, "Novelty Image");

    timer.stop();
    cout << "AL time for writing the raw novelty image: " << timer.getLastAbsolute() << endl;
  }

}

inline void SemSegNovelty::computeClassificationResults( const NICE::MultiChannelImageT<double> & feats, 
                                                   NICE::ImageT<int> & segresult,
                                                   NICE::MultiChannelImageT<double> & probabilities,
                                                   const int & xsize,
                                                   const int & ysize,
                                                   const int & featdim
                                                       )
{
  std::cerr << "featdim: " << featdim << std::endl;
  
  if ( classifier != NULL )
  {  
  
    #pragma omp parallel for
    for ( int y = 0; y < ysize; y += testWSize )
    {
      Example example;
      example.vec = NULL;
      example.svec = new SparseVector ( featdim );
      for ( int x = 0; x < xsize; x += testWSize)
      {
        for ( int f = 0; f < featdim; f++ )
        {
          double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
          if ( val > 1e-10 )
            ( *example.svec ) [f] = val;
        }
        example.svec->normalize();

        ClassificationResult cr = classifier->classify ( example );

        int xs = std::max(0, x - testWSize/2);
        int xe = std::min(xsize - 1, x + testWSize/2);
        int ys = std::max(0, y - testWSize/2);
        int ye = std::min(ysize - 1, y + testWSize/2);
        for (int yl = ys; yl <= ye; yl++)
        {
          for (int xl = xs; xl <= xe; xl++)
          {
            for ( int j = 0 ; j < cr.scores.size(); j++ )
            {
              probabilities ( xl, yl, j ) = cr.scores[j];
            }
            segresult ( xl, yl ) = cr.classno;
          }
        }
        
        example.svec->clear();
      }
      delete example.svec;
      example.svec = NULL;
    }
  }
  else //vclassifier
  {
    std::cerr << "compute classification results with vclassifier" << std::endl;
    #pragma omp parallel for
    for ( int y = 0; y < ysize; y += testWSize )
    {
      for ( int x = 0; x < xsize; x += testWSize)
      {
        NICE::Vector v(featdim);
        for ( int f = 0; f < featdim; f++ )
        {
          double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
          v[f] = val;
        }
        v.normalizeL1();

        ClassificationResult cr = vclassifier->classify ( v );

        int xs = std::max(0, x - testWSize/2);
        int xe = std::min(xsize - 1, x + testWSize/2);
        int ys = std::max(0, y - testWSize/2);
        int ye = std::min(ysize - 1, y + testWSize/2);
        for (int yl = ys; yl <= ye; yl++)
        {
          for (int xl = xs; xl <= xe; xl++)
          {
            for ( int j = 0 ; j < cr.scores.size(); j++ )
            {
              probabilities ( xl, yl, j ) = cr.scores[j];
            }
            segresult ( xl, yl ) = cr.classno;
          }
        }
      }
    }    

  }
}

// compute novelty images depending on the strategy chosen

void SemSegNovelty::computeNoveltyByRandom(         NICE::FloatImage & noveltyImage, 
                                              const NICE::MultiChannelImageT<double> & feats,  
                                                    NICE::ImageT<int> & segresult,
                                                    NICE::MultiChannelImageT<double> & probabilities,
                                             const int & xsize, const int & ysize, const int & featdim  )
{
#pragma omp parallel for
  for ( int y = 0; y < ysize; y += testWSize )
  {
    Example example;
    example.vec = NULL;
    example.svec = new SparseVector ( featdim );
    for ( int x = 0; x < xsize; x += testWSize)
    {
      for ( int f = 0; f < featdim; f++ )
      {
        double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
        if ( val > 1e-10 )
          ( *example.svec ) [f] = val;
      }
      example.svec->normalize();

      ClassificationResult cr = classifier->classify ( example );
  
      int xs = std::max(0, x - testWSize/2);
      int xe = std::min(xsize - 1, x + testWSize/2);
      int ys = std::max(0, y - testWSize/2);
      int ye = std::min(ysize - 1, y + testWSize/2);
      
      double randVal = randDouble();

      for (int yl = ys; yl <= ye; yl++)
      {
        for (int xl = xs; xl <= xe; xl++)
        {
          for ( int j = 0 ; j < cr.scores.size(); j++ )
          {
            probabilities ( xl, yl, j ) = cr.scores[j];
          }
          segresult ( xl, yl ) = cr.classno;          
          noveltyImage ( xl, yl ) = randVal; 
        }
      }     
    }
  }  
}


void SemSegNovelty::computeNoveltyByVariance(       NICE::FloatImage & noveltyImage, 
                                              const NICE::MultiChannelImageT<double> & feats,  
                                                    NICE::ImageT<int> & segresult,
                                                    NICE::MultiChannelImageT<double> & probabilities,
                                             const int & xsize, const int & ysize, const int & featdim )
{
#pragma omp parallel for
  for ( int y = 0; y < ysize; y += testWSize )
  {
    Example example;
    example.vec = NULL;
    example.svec = new SparseVector ( featdim );
    for ( int x = 0; x < xsize; x += testWSize)
    {
      for ( int f = 0; f < featdim; f++ )
      {
        double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
        if ( val > 1e-10 )
          ( *example.svec ) [f] = val;
      }
      example.svec->normalize();

      ClassificationResult cr = classifier->classify ( example );
      
      int xs = std::max(0, x - testWSize/2);
      int xe = std::min(xsize - 1, x + testWSize/2);
      int ys = std::max(0, y - testWSize/2);
      int ye = std::min(ysize - 1, y + testWSize/2);
      for (int yl = ys; yl <= ye; yl++)
      {
        for (int xl = xs; xl <= xe; xl++)
        {
          for ( int j = 0 ; j < cr.scores.size(); j++ )
          {
            probabilities ( xl, yl, j ) = cr.scores[j];
          }
          segresult ( xl, yl ) = cr.classno;          
          noveltyImage ( xl, yl ) = cr.uncertainty; 
        }
      }          
      
      example.svec->clear();
    }
    delete example.svec;
    example.svec = NULL;
  }  
}

void SemSegNovelty::computeNoveltyByGPUncertainty(  NICE::FloatImage & noveltyImage, 
                                              const NICE::MultiChannelImageT<double> & feats,  
                                                    NICE::ImageT<int> & segresult,
                                                    NICE::MultiChannelImageT<double> & probabilities,
                                             const int & xsize, const int & ysize, const int & featdim )
{
  
  double gpNoise =  0.01;
  //TODO getMethod for GPHIK
  //conf->gD("GPHIK", "noise", 0.01);
  
  
#pragma omp parallel for
  for ( int y = 0; y < ysize; y += testWSize )
  {
    Example example;
    example.vec = NULL;
    example.svec = new SparseVector ( featdim );
    for ( int x = 0; x < xsize; x += testWSize)
    {
      for ( int f = 0; f < featdim; f++ )
      {
        double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
        if ( val > 1e-10 )
          ( *example.svec ) [f] = val;
      }
      example.svec->normalize();

      ClassificationResult cr = classifier->classify ( example );
      
      double maxMeanAbs ( 0.0 ); 
      
      for ( int j = 0 ; j < cr.scores.size(); j++ )
      {   
        if ( forbidden_classesTrain.find ( j ) != forbidden_classesTrain.end() )
        {
          continue;
        }
        //check for larger abs mean
        if (abs(cr.scores[j]) > maxMeanAbs)
        {
          maxMeanAbs = abs(cr.scores[j]);
        }
        
      }

      double firstTerm (1.0 / sqrt(cr.uncertainty+gpNoise));
      
      //compute the heuristic GP-UNCERTAINTY, as proposed by Kapoor et al. in IJCV 2010
      // GP-UNCERTAINTY : |mean| / sqrt(var^2 + gpnoise^2)
      double gpUncertaintyVal = maxMeanAbs*firstTerm; //firstTerm = 1.0 / sqrt(r.uncertainty+gpNoise))

      int xs = std::max(0, x - testWSize/2);
      int xe = std::min(xsize - 1, x + testWSize/2);
      int ys = std::max(0, y - testWSize/2);
      int ye = std::min(ysize - 1, y + testWSize/2);
      for (int yl = ys; yl <= ye; yl++)
      {
        for (int xl = xs; xl <= xe; xl++)
        {         
          for ( int j = 0 ; j < cr.scores.size(); j++ )
          {
            probabilities ( xl, yl, j ) = cr.scores[j];
          }
          segresult ( xl, yl ) = cr.classno;          
          noveltyImage ( xl, yl ) = gpUncertaintyVal;  
        }
      }   
      
      example.svec->clear();
    }
    delete example.svec;
    example.svec = NULL;
  }  
}

void SemSegNovelty::computeNoveltyByGPMean(  NICE::FloatImage & noveltyImage, 
                                              const NICE::MultiChannelImageT<double> & feats,  
                                                    NICE::ImageT<int> & segresult,
                                                    NICE::MultiChannelImageT<double> & probabilities,
                                             const int & xsize, const int & ysize, const int & featdim )
{
  double gpNoise =  0.01;
  //TODO getMethod for GPHIK
  //conf->gD("GPHIK", "noise", 0.01);
  
#pragma omp parallel for
  for ( int y = 0; y < ysize; y += testWSize )
  {
    Example example;
    example.vec = NULL;
    example.svec = new SparseVector ( featdim );
    for ( int x = 0; x < xsize; x += testWSize)
    {
      for ( int f = 0; f < featdim; f++ )
      {
        double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
        if ( val > 1e-10 )
          ( *example.svec ) [f] = val;
      }
      example.svec->normalize();

      ClassificationResult cr = classifier->classify ( example );

      double minMeanAbs ( numeric_limits<double>::max() );
      
      for ( int j = 0 ; j < probabilities.channels(); j++ )
      {
        if ( forbidden_classesTrain.find ( j ) != forbidden_classesTrain.end() )
        {
          continue;
        }

        //check whether we found a class with higher smaller abs mean than the current minimum
        if (abs( cr.scores[j] ) < minMeanAbs)
        {
          minMeanAbs = abs(cr.scores[j]); 
        }
      }

      // compute results when we take the lowest mean value of all classes
      double gpMeanVal = minMeanAbs;
  
      int xs = std::max(0, x - testWSize/2);
      int xe = std::min(xsize - 1, x + testWSize/2);
      int ys = std::max(0, y - testWSize/2);
      int ye = std::min(ysize - 1, y + testWSize/2);
      for (int yl = ys; yl <= ye; yl++)
      {
        for (int xl = xs; xl <= xe; xl++)
        {
          for ( int j = 0 ; j < cr.scores.size(); j++ )
          {
            probabilities ( xl, yl, j ) = cr.scores[j];
          }
          segresult ( xl, yl ) = cr.classno;          
          noveltyImage ( xl, yl ) = gpMeanVal; 
        }
      }     
    }
  }  
}

void SemSegNovelty::computeNoveltyByGPMeanRatio(  NICE::FloatImage & noveltyImage, 
                                              const NICE::MultiChannelImageT<double> & feats,  
                                                    NICE::ImageT<int> & segresult,
                                                    NICE::MultiChannelImageT<double> & probabilities,
                                             const int & xsize, const int & ysize, const int & featdim )
{
  double gpNoise =  0.01;
  //TODO getMethod for GPHIK
  //conf->gD("GPHIK", "noise", 0.01);
  
#pragma omp parallel for
  for ( int y = 0; y < ysize; y += testWSize )
  {
    Example example;
    example.vec = NULL;
    example.svec = new SparseVector ( featdim );
    for ( int x = 0; x < xsize; x += testWSize)
    {
      for ( int f = 0; f < featdim; f++ )
      {
        double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
        if ( val > 1e-10 )
          ( *example.svec ) [f] = val;
      }
      example.svec->normalize();

      ClassificationResult cr = classifier->classify ( example );

      double maxMean ( -numeric_limits<double>::max() );
      double sndMaxMean ( -numeric_limits<double>::max() );     
      
      for ( int j = 0 ; j < cr.scores.size(); j++ )
      {
        if ( forbidden_classesTrain.find ( j ) != forbidden_classesTrain.end() )
        {
          continue;
        }
        
        //check for larger mean without abs as well
        if (cr.scores[j] > maxMean)
        {
          sndMaxMean = maxMean;
          maxMean = cr.scores[j];
        }
        // and also for the second highest mean of all classes
        else if (cr.scores[j] > sndMaxMean)
        {
          sndMaxMean = cr.scores[j];
        }          
      }
      
      //look at the difference in the absolut mean values for the most plausible class
      // and the second most plausible class
      double gpMeanRatioVal= maxMean - sndMaxMean;


      int xs = std::max(0, x - testWSize/2);
      int xe = std::min(xsize - 1, x + testWSize/2);
      int ys = std::max(0, y - testWSize/2);
      int ye = std::min(ysize - 1, y + testWSize/2);
      for (int yl = ys; yl <= ye; yl++)
      {
        for (int xl = xs; xl <= xe; xl++)
        {
          for ( int j = 0 ; j < cr.scores.size(); j++ )
          {
            probabilities ( xl, yl, j ) = cr.scores[j];
          }
          segresult ( xl, yl ) = cr.classno;          
          noveltyImage ( xl, yl ) = gpMeanRatioVal;
        }
      }    
      
      example.svec->clear();
    }
    delete example.svec;
    example.svec = NULL;
  }  
}

void SemSegNovelty::computeNoveltyByGPWeightAll(  NICE::FloatImage & noveltyImage, 
                                              const NICE::MultiChannelImageT<double> & feats,  
                                                    NICE::ImageT<int> & segresult,
                                                    NICE::MultiChannelImageT<double> & probabilities,
                                             const int & xsize, const int & ysize, const int & featdim )
{
  double gpNoise =  0.01;
  //TODO getMethod for GPHIK
  //conf->gD("GPHIK", "noise", 0.01);
  
#pragma omp parallel for
  for ( int y = 0; y < ysize; y += testWSize )
  {
    Example example;
    example.vec = NULL;
    example.svec = new SparseVector ( featdim );
    for ( int x = 0; x < xsize; x += testWSize)
    {
      for ( int f = 0; f < featdim; f++ )
      {
        double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
        if ( val > 1e-10 )
          ( *example.svec ) [f] = val;
      }
      example.svec->normalize();

      ClassificationResult cr = classifier->classify ( example );
      
      double firstTerm (1.0 / sqrt(cr.uncertainty+gpNoise));
      
      double gpWeightAllVal ( 0.0 );

      if ( numberOfClasses > 2)
      {
        //compute the weight in the alpha-vector for every sample after assuming it to be 
        // added to the training set.
        // Thereby, we measure its "importance" for the current model
        // 
        //double firstTerm is already computed
        //
        //the second term is only needed when computing impacts
        //double secondTerm; //this is the nasty guy :/
        
        //--- compute the third term
        // this is the difference between predicted label and GT label 
        std::vector<double> diffToPositive; diffToPositive.clear();
        std::vector<double> diffToNegative; diffToNegative.clear();
        double diffToNegativeSum(0.0);
        
        for ( int j = 0 ; j < cr.scores.size(); j++ )
        {
          if ( forbidden_classesTrain.find ( j ) != forbidden_classesTrain.end() )
          {
            continue;
          }          
          
          // look at the difference to plus 1          
          diffToPositive.push_back(abs(cr.scores[j] - 1));
          // look at the difference to -1          
          diffToNegative.push_back(abs(cr.scores[j] + 1));
          //sum up the difference to -1
          diffToNegativeSum += abs(cr.scores[j] - 1);
        }

        //let's subtract for every class its diffToNegative from the sum, add its diffToPositive,
        //and use this as the third term for this specific class.
        //the final value is obtained by minimizing over all classes
        //
        // originally, we minimize over all classes after building the final score
        // however, the first and the second term do not depend on the choice of
        // y*, therefore we minimize here already
        double thirdTerm (numeric_limits<double>::max()) ;
        for(uint tmpCnt = 0; tmpCnt < diffToPositive.size(); tmpCnt++)
        {
          double tmpVal ( diffToPositive[tmpCnt] + (diffToNegativeSum-diffToNegative[tmpCnt])   );
          if (tmpVal < thirdTerm)
            thirdTerm = tmpVal;
        }
        gpWeightAllVal = thirdTerm*firstTerm;        
      }
      else //binary scenario
      {
        gpWeightAllVal = std::min( abs(cr.scores[*classesInUse.begin()]+1), abs(cr.scores[*classesInUse.begin()]-1) );
        gpWeightAllVal *= firstTerm;
      }

      int xs = std::max(0, x - testWSize/2);
      int xe = std::min(xsize - 1, x + testWSize/2);
      int ys = std::max(0, y - testWSize/2);
      int ye = std::min(ysize - 1, y + testWSize/2);
      for (int yl = ys; yl <= ye; yl++)
      {
        for (int xl = xs; xl <= xe; xl++)
        {
          for ( int j = 0 ; j < cr.scores.size(); j++ )
          {
            probabilities ( xl, yl, j ) = cr.scores[j];
          }
          segresult ( xl, yl ) = cr.classno;          
          noveltyImage ( xl, yl ) = gpWeightAllVal;
        }
      }
   
      
      example.svec->clear();
    }
    delete example.svec;
    example.svec = NULL;
  }  
}

void SemSegNovelty::computeNoveltyByGPWeightRatio(  NICE::FloatImage & noveltyImage, 
                                              const NICE::MultiChannelImageT<double> & feats,  
                                                    NICE::ImageT<int> & segresult,
                                                    NICE::MultiChannelImageT<double> & probabilities,
                                             const int & xsize, const int & ysize, const int & featdim )
{
  double gpNoise =  0.01;
  //TODO getMethod for GPHIK
  //conf->gD("GPHIK", "noise", 0.01);
  
#pragma omp parallel for
  for ( int y = 0; y < ysize; y += testWSize )
  {
    Example example;
    example.vec = NULL;
    example.svec = new SparseVector ( featdim );
    for ( int x = 0; x < xsize; x += testWSize)
    {
      for ( int f = 0; f < featdim; f++ )
      {
        double val = feats.getIntegralValue ( x - whs, y - whs, x + whs, y + whs, f );
        if ( val > 1e-10 )
          ( *example.svec ) [f] = val;
      }
      example.svec->normalize();

      ClassificationResult cr = classifier->classify ( example );
 

       double firstTerm (1.0 / sqrt(cr.uncertainty+gpNoise));

       double gpWeightRatioVal ( 0.0 );

       if ( numberOfClasses > 2)
       {
        //compute the weight in the alpha-vector for every sample after assuming it to be 
        // added to the training set.
        // Thereby, we measure its "importance" for the current model
        // 
        //double firstTerm is already computed
        //
        //the second term is only needed when computing impacts
        //double secondTerm; //this is the nasty guy :/
        
        //--- compute the third term
        // this is the difference between predicted label and GT label 
        std::vector<double> diffToPositive; diffToPositive.clear();
        std::vector<double> diffToNegative; diffToNegative.clear();
        double diffToNegativeSum(0.0);
        
        for ( int j = 0 ; j < cr.scores.size(); j++ )
        {
          if ( forbidden_classesTrain.find ( j ) != forbidden_classesTrain.end() )
          {
            continue;
          }          
          
          // look at the difference to plus 1          
          diffToPositive.push_back(abs(cr.scores[j] - 1));
        }

        //let's subtract for every class its diffToNegative from the sum, add its diffToPositive,
        //and use this as the third term for this specific class.
        //the final value is obtained by minimizing over all classes
        //
        // originally, we minimize over all classes after building the final score
        // however, the first and the second term do not depend on the choice of
        // y*, therefore we minimize here already   
        
        //now look on the ratio of the resulting weights for the most plausible
        // against the second most plausible class
        double thirdTermMostPlausible ( 0.0 ) ;
        double thirdTermSecondMostPlausible ( 0.0 ) ;
        for(uint tmpCnt = 0; tmpCnt < diffToPositive.size(); tmpCnt++)
        {
          if (diffToPositive[tmpCnt] > thirdTermMostPlausible)
          {
            thirdTermSecondMostPlausible = thirdTermMostPlausible;
            thirdTermMostPlausible = diffToPositive[tmpCnt];
          }
          else if (diffToPositive[tmpCnt] > thirdTermSecondMostPlausible)
          {
            thirdTermSecondMostPlausible = diffToPositive[tmpCnt];
          }
        }
        //compute the resulting score
        gpWeightRatioVal = (thirdTermMostPlausible - thirdTermSecondMostPlausible)*firstTerm;      

        //finally, look for this feature how it would affect to whole model (summarized by weight-vector alpha), if we would 
        //use it as an additional training example
        //TODO this would be REALLY computational demanding. Do we really want to do this?
  //         gpImpactAll[s] ( pce[i].second.x, pce[i].second.y ) = thirdTerm*firstTerm*secondTerm;
  //         gpImpactRatio[s] ( pce[i].second.x, pce[i].second.y ) = (thirdTermMostPlausible - thirdTermSecondMostPlausible)*firstTerm*secondTerm;
       }
       else //binary scenario
       {
         gpWeightRatioVal = std::min( abs(cr.scores[*classesInUse.begin()]+1), abs(cr.scores[*classesInUse.begin()]-1) );
         gpWeightRatioVal *= firstTerm;
       }

      int xs = std::max(0, x - testWSize/2);
      int xe = std::min(xsize - 1, x + testWSize/2);
      int ys = std::max(0, y - testWSize/2);
      int ye = std::min(ysize - 1, y + testWSize/2);
      for (int yl = ys; yl <= ye; yl++)
      {
        for (int xl = xs; xl <= xe; xl++)
        {
          for ( int j = 0 ; j < cr.scores.size(); j++ )
          {
            probabilities ( xl, yl, j ) = cr.scores[j];
          }
          segresult ( xl, yl ) = cr.classno;          
          noveltyImage ( xl, yl ) = gpWeightRatioVal;  
        }
      }
       
      example.svec->clear();
    }
    delete example.svec;
    example.svec = NULL;
  }  
}


void SemSegNovelty::addNewExample(const NICE::Vector& v_newExample, const int & newClassNo)
{
  //accept the new class as valid information
  if ( forbidden_classesTrain.find ( newClassNo ) != forbidden_classesTrain.end() )
  {
    forbidden_classesTrain.erase(newClassNo);
    numberOfClasses++;
  }
  if ( classesInUse.find ( newClassNo ) == classesInUse.end() )
  {
    classesInUse.insert( newClassNo );
  }    
  
  
  //then add it to the classifier used
  if ( classifier != NULL )
  { 
    if (this->classifierString.compare("GPHIKClassifier") == 0)    
    {
      Example newExample;
      SparseVector svec ( v_newExample );
      newExample.svec = &svec;
      static_cast<FPCGPHIK*>(classifier)->addExample ( newExample, newClassNo );
    } 
  }
  else //vclassifier
  {
    if (this->classifierString.compare("nn") == 0)    
    {
      vclassifier->teach ( newClassNo, v_newExample );
    }
  }
}

void SemSegNovelty::addNovelExamples()
{

  Timer timer;
  
  //show the image that contains the most novel region
  if (b_visualizeALimages)
    showImage(maskedImg, "Most novel region");  
  
  timer.start();
  

    std::stringstream out;
    std::vector< std::string > list2;
    StringTools::split ( Globals::getCurrentImgFN (), '/', list2 );
    out << resultdir << "/" << list2.back();     
    
    maskedImg.writePPM ( out.str() + "_run_" +  NICE::intToString(this->iterationCountSuffix) + "_" + noveltyMethodString+ "_query.ppm" );

  
  timer.stop();
  std::cerr << "AL time for writing queried image: " << timer.getLast() << std::endl;

  timer.start();
  
  //check which classes will be added using the features from the novel region
  std::set<int> newClassNumbers;
  newClassNumbers.clear(); //just to be sure  
  for ( uint i = 0 ; i < newTrainExamples.size() ; i++ )
  {
    if (newClassNumbers.find(newTrainExamples[i].first /* classNumber*/) == newClassNumbers.end() )
    {
      newClassNumbers.insert(newTrainExamples[i].first );
    }
  }

  //accept the new classes as valid information
  for (std::set<int>::const_iterator clNoIt = newClassNumbers.begin(); clNoIt != newClassNumbers.end(); clNoIt++)
  {
    if ( forbidden_classesTrain.find ( *clNoIt ) != forbidden_classesTrain.end() )
    {
      forbidden_classesTrain.erase(*clNoIt);
      numberOfClasses++;
    }
    if ( classesInUse.find ( *clNoIt ) == classesInUse.end() )
    {
      classesInUse.insert( *clNoIt );
    }
  }
  
  timer.stop();
  std::cerr << "AL time for accepting possible new classes: " << timer.getLast() << std::endl;
  
  timer.start();
  //then add the new features to the classifier used
  if ( classifier != NULL )
  { 
    if (this->classifierString.compare("GPHIKClassifier") == 0)    
    {
      classifier->addMultipleExamples ( this->newTrainExamples );
    }    
  }
  else //vclassifier
  {
    //TODO
  }
  
  timer.stop();
  std::cerr << "AL time for actually updating the classifier: " << timer.getLast() << std::endl;
  
  std::cerr << "the current region to query is: " << currentRegionToQuery.first << " -- " << currentRegionToQuery.second << std::endl;
  
  //did we already query a region of this image?
  if ( queriedRegions.find( currentRegionToQuery.first ) != queriedRegions.end() )
  {
    queriedRegions[ currentRegionToQuery.first ].insert(currentRegionToQuery.second);
  }
  else
  {
    std::set<int> tmpSet; tmpSet.insert(currentRegionToQuery.second);
    queriedRegions.insert(std::pair<std::string,std::set<int> > (currentRegionToQuery.first, tmpSet ) );
  }  
  
  std::cerr << "Write already queried regions: " << std::endl;
  for (std::map<std::string,std::set<int> >::const_iterator it = queriedRegions.begin(); it != queriedRegions.end(); it++)
  {
    std::cerr << "image: " << it->first << " --   ";
    for (std::set<int>::const_iterator itReg = it->second.begin(); itReg != it->second.end(); itReg++)
    {
      std::cerr << *itReg << " ";
    } 
    std::cerr << std::endl;
  }
  
  //clear the latest results, since one iteration is over
  globalMaxUncert = -numeric_limits<double>::max();
  if (!mostNoveltyWithMaxScores)
    globalMaxUncert = numeric_limits<double>::max();
}

const Examples * SemSegNovelty::getNovelExamples() const
{
  return &(this->newTrainExamples);
}

///////////////////// INTERFACE PERSISTENT /////////////////////
// interface specific methods for store and restore
///////////////////// INTERFACE PERSISTENT ///////////////////// 

void SemSegNovelty::restore ( std::istream & is, int format )
{
  //delete everything we knew so far...
  this->clear();
  
  bool b_restoreVerbose ( false );
#ifdef B_RESTOREVERBOSE
  b_restoreVerbose = true;
#endif  
  
  if ( is.good() )
  {
    if ( b_restoreVerbose ) 
      std::cerr << " restore SemSegNovelty" << std::endl;
    
    std::string tmp;
    is >> tmp; //class name 
    
    if ( ! this->isStartTag( tmp, "SemSegNovelty" ) )
    {
      std::cerr << " WARNING - attempt to restore SemSegNovelty, but start flag " << tmp << " does not match! Aborting... " << std::endl;
      throw;
    }   
    
    if (classifier != NULL)
    {
      delete classifier;
      classifier = NULL;
    }    
    
    is.precision (numeric_limits<double>::digits10 + 1);
    
    bool b_endOfBlock ( false ) ;
    
    while ( !b_endOfBlock )
    {
      is >> tmp; // start of block 
      
      if ( this->isEndTag( tmp, "SemSegNovelty" ) )
      {
        b_endOfBlock = true;
        continue;
      }      
      
      tmp = this->removeStartTag ( tmp );
      
      if ( b_restoreVerbose )
        std::cerr << " currently restore section " << tmp << " in SemSegNovelty" << std::endl;
      
      
      ///////////////////////////////
      //     FEATURE EXTRACTION    //
      ///////////////////////////////
      if ( tmp.compare("featExtract") == 0 )
      { 
        featExtract->restore(is, format);
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }      
      else if ( tmp.compare("trainWSize") == 0 )
      { 
        is >> trainWSize;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("whs") == 0 )
      { 
        is >> whs;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("testWSize") == 0 )
      { 
        is >> testWSize;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }  
      ///////////////////////////////
      //     NOVELTY COMPUTATION   //
      ///////////////////////////////
      else if ( tmp.compare("noveltyMethod") == 0 )
      {
        unsigned int ui_noveltyMethod;
        is >> ui_noveltyMethod;        
        this->noveltyMethod = static_cast<NoveltyMethod> ( ui_noveltyMethod );
	
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("noveltyMethodString") == 0 )
      {
        is >> noveltyMethodString;        
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("globalMaxUncert") == 0 )
      {
        is >> globalMaxUncert;        
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("mostNoveltyWithMaxScores") == 0 )
      { 
        is >> mostNoveltyWithMaxScores;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("findMaximumUncert") == 0 )
      { 
        is >> findMaximumUncert;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      //TODO maskedImg
      else if ( tmp.compare("b_visualizeALimages") == 0 )
      { 
        is >> b_visualizeALimages;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      ///////////////////////////////
      //     CLASSIFICATION STUFF  //
      /////////////////////////////// 
      else if ( tmp.compare("classifier") == 0 )
      {
        std::string isNull;
        is >> isNull;
        
        // check whether we originally used a classifier
        if ( isNull.compare( "NULL" ) == 0 )
        {
          if ( classifier != NULL )
            delete classifier;
          classifier = NULL;
        }
        else
        {
          if ( classifier == NULL )
            classifier = new OBJREC::FPCGPHIK();

          classifier->restore(is, format);              
        }
          
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("vclassifier") == 0 )
      {
        std::string isNull;
        is >> isNull;
        
        // check whether we originally used a vclassifier
        if ( isNull.compare( "NULL" ) == 0 )
        {
          if ( vclassifier != NULL )
            delete vclassifier;
          vclassifier = NULL;
        }
        else
        {
          fthrow ( NICE::Exception, "Restoring of VecClassifiers is not implemented yet!" );
/*          if ( vclassifier == NULL )
            vclassifier = new OBJREC::VecClassifier();

          vclassifier->restore(is, format);  */          
        }
          
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("forbidden_classesTrain") == 0 )
      {
        is >> tmp; // size
        int forbClTrainSize ( 0 );
        is >> forbClTrainSize;

        forbidden_classesTrain.clear();
        
        if ( b_restoreVerbose ) 
          std::cerr << "restore forbidden_classesTrain with size: " << forbClTrainSize << std::endl;

        if ( forbClTrainSize > 0 )
        {
          if ( b_restoreVerbose ) 
            std::cerr << " restore forbidden_classesTrain" << std::endl;
          
          for (int i = 0; i < forbClTrainSize; i++)
          {
            int classNo;
            is >> classNo;        
            forbidden_classesTrain.insert ( classNo );
          }
        } 
        else
        {
          if ( b_restoreVerbose ) 
            std::cerr << " skip restoring forbidden_classesTrain" << std::endl;
        }
        
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("forbidden_classesActiveLearning") == 0 )
      {
        is >> tmp; // size
        int forbClALSize ( 0 );
        is >> forbClALSize;

        forbidden_classesActiveLearning.clear();
        
        if ( b_restoreVerbose ) 
          std::cerr << "restore forbidden_classesActiveLearning with size: " << forbClALSize << std::endl;

        if ( forbClALSize > 0 )
        {
          if ( b_restoreVerbose ) 
            std::cerr << " restore forbidden_classesActiveLearning" << std::endl;
          
          for (int i = 0; i < forbClALSize; i++)
          {
            int classNo;
            is >> classNo;        
            forbidden_classesActiveLearning.insert ( classNo );
          }
        } 
        else
        {
          if ( b_restoreVerbose ) 
            std::cerr << " skip restoring forbidden_classesActiveLearning" << std::endl;
        }
        
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("classesInUse") == 0 )
      {
        is >> tmp; // size
        int clInUseSize ( 0 );
        is >> clInUseSize;

        classesInUse.clear();
        
        if ( b_restoreVerbose ) 
          std::cerr << "restore classesInUse with size: " << clInUseSize << std::endl;

        if ( clInUseSize > 0 )
        {
          if ( b_restoreVerbose ) 
            std::cerr << " restore classesInUse" << std::endl;
          
          for (int i = 0; i < clInUseSize; i++)
          {
            int classNo;
            is >> classNo;        
            classesInUse.insert ( classNo );
          }
        } 
        else
        {
          if ( b_restoreVerbose ) 
            std::cerr << " skip restoring classesInUse" << std::endl;
        }
        
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("numberOfClasses") == 0 )
      { 
        is >> numberOfClasses;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("read_classifier") == 0 )
      { 
        is >> read_classifier;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("save_classifier") == 0 )
      { 
        is >> save_classifier;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      } 
      else if ( tmp.compare("cache") == 0 )
      { 
        is >> cache;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      else if ( tmp.compare("resultdir") == 0 )
      { 
        is >> resultdir;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      //TODO newTrainExamples
      ///////////////////////////////
      //     SEGMENTATION STUFF    //
      ///////////////////////////////
      //TODO regionSeg
      else if ( tmp.compare("s_rsMethode") == 0 )
      { 
        is >> this->s_rsMethode;
	// theoretically, we should properly store and restore the regionSeg object. However,  its parent class does not provide
	// a Persistent interface yet. Hence, we perform this tiny workaround which works, since regionSeg is not changed over time...
	// only be aware of parameters originally set via config...
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }
      //NOTE regionSeg seems really important to keep track off
      else if ( tmp.compare("reuseSegmentation") == 0 )
      { 
        is >> reuseSegmentation;
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }       
      else if ( tmp.compare("queriedRegions") == 0 )
      {
        is >> tmp; // size
        int queriedRegionsSize ( 0 );
        is >> queriedRegionsSize;
        queriedRegions.clear();

        if ( b_restoreVerbose ) 
          std::cerr << "restore queriedRegions with size: " << queriedRegionsSize << std::endl;
        for ( int i = 0; i < queriedRegionsSize; i++ )
        {
	  // restore key
          std::string key;
          is >> key;
	  
	  // restore values -- inner loop over sets
	  is >> tmp; // size
	  int regionsOfImgSize ( 0 );
	  is >> regionsOfImgSize;

	  std::set< int > regionsOfImg;
	  regionsOfImg.clear();
	    
	    for (int i = 0; i < regionsOfImgSize; i++)
	    {
	      int idxRegion;
	      is >> idxRegion;        
	      regionsOfImg.insert ( idxRegion );
	    }
          queriedRegions.insert ( std::pair<std::string, std::set< int > > ( key, regionsOfImg ) );
        }
        
        is >> tmp; // end of block 
        tmp = this->removeEndTag ( tmp );
      }            
      //
      //TODO currentRegionToQuery
      //
      ///////////////////////////////
      //       PARENT OBJECT       //
      ///////////////////////////////
      else if ( tmp.compare("SemSegNovelty--Parent") == 0 )
      {
        // restore parent object
        SemanticSegmentation::restore(is);
      }      
      else
      {
      std::cerr << "WARNING -- unexpected SemSegNovelty object -- " << tmp << " -- for restoration... aborting" << std::endl;
      throw;
      }
      
      // INSTANTIATE (YET) NON-RESTORABLE OBJECTS
      //TODO destructor of regionSeg is non-virtual so far - change this accordingly!
      if ( this->regionSeg != NULL ) 
	delete this->regionSeg;
      
      if( this->s_rsMethode == "none" )
      {
	this->regionSeg = NULL;
      }
      else
      {
	//NOTE using an empty config file might not be save...
	NICE::Config tmpConfEmpty;
	RegionSegmentationMethod *tmpRegionSeg = GenericRegionSegmentationMethodSelection::selectRegionSegmentationMethod( &tmpConfEmpty, this->s_rsMethode );    
	if ( reuseSegmentation )
	  this->regionSeg = new RSCache ( &tmpConfEmpty, tmpRegionSeg );
	else
	  this->regionSeg = tmpRegionSeg;
      }      
      
      // done restoration
    }
  }
  else
  {
    std::cerr << "SemSegNovelty::restore -- InStream not initialized - restoring not possible!" << std::endl;
    throw;
  }
  
  
 
}

void SemSegNovelty::store ( std::ostream & os, int format ) const
{ 
  if (os.good())
  {
    // show starting point
    os << this->createStartTag( "SemSegNovelty" ) << std::endl;    
    
    ///////////////////////////////
    //     FEATURE EXTRACTION    //
    ///////////////////////////////
    os << this->createStartTag( "featExtract" ) << std::endl;
    featExtract->store ( os );
    os << this->createStartTag( "featExtract" ) << std::endl; 
    
    os << this->createStartTag( "trainWsize" ) << std::endl;
    os << this->trainWSize << std::endl;
    os << this->createStartTag( "trainWsize" ) << std::endl;  
    
    os << this->createStartTag( "whs" ) << std::endl;
    os << this->whs << std::endl;
    os << this->createStartTag( "whs" ) << std::endl;
    
    os << this->createStartTag( "testWSize" ) << std::endl;
    os << this->testWSize << std::endl;
    os << this->createStartTag( "testWSize" ) << std::endl;      
    
    ///////////////////////////////
    //     NOVELTY COMPUTATION   //
    ///////////////////////////////
    
    os << this->createStartTag( "noveltyMethod" ) << std::endl;
    os << this->noveltyMethod << std::endl;
    os << this->createStartTag( "noveltyMethod" ) << std::endl;

    os << this->createStartTag( "noveltyMethodString" ) << std::endl;
    os << this->noveltyMethodString << std::endl;
    os << this->createStartTag( "noveltyMethodString" ) << std::endl;
    
    os << this->createStartTag( "globalMaxUncert" ) << std::endl;
    os << this->globalMaxUncert << std::endl;
    os << this->createStartTag( "globalMaxUncert" ) << std::endl;

    os << this->createStartTag( "mostNoveltyWithMaxScores" ) << std::endl;
    os << this->mostNoveltyWithMaxScores << std::endl;
    os << this->createStartTag( "mostNoveltyWithMaxScores" ) << std::endl;

    os << this->createStartTag( "findMaximumUncert" ) << std::endl;
    os << this->findMaximumUncert << std::endl;
    os << this->createStartTag( "findMaximumUncert" ) << std::endl;    

    //TODO maskedImg
    
    os << this->createStartTag( "b_visualizeALimages" ) << std::endl;
    os << this->b_visualizeALimages << std::endl;
    os << this->createStartTag( "b_visualizeALimages" ) << std::endl;  

    
    ///////////////////////////////
    //     CLASSIFICATION STUFF  //
    ///////////////////////////////
    
    os << this->createStartTag( "classifierString" ) << std::endl;
    os << this->classifierString << std::endl;
    os << this->createStartTag( "classifierString" ) << std::endl; 

    os << this->createStartTag( "classifier" ) << std::endl;
    if ( this->classifier != NULL )
    {
      os << "NOTNULL" << std::endl;
      classifier->store ( os, format );
    }
    else
    {
     os << "NULL" << std::endl;
    }
    os << this->createEndTag( "classifier" ) << std::endl;
    
    //
    
    os << this->createStartTag( "vclassifier" ) << std::endl;
    if ( this->classifier != NULL )
    {
      os << "NOTNULL" << std::endl;
      vclassifier->store ( os, format );
    }
    else
    {
     os << "NULL" << std::endl;
    }
    os << this->createEndTag( "vclassifier" ) << std::endl;
    
    
    os << this->createStartTag( "forbidden_classesTrain" ) << std::endl;
    os << "size: " << forbidden_classesTrain.size() << std::endl;

    for ( std::set< int >::const_iterator itForbClassTrain = forbidden_classesTrain.begin();
          itForbClassTrain != forbidden_classesTrain.end();
          itForbClassTrain++
        )
    {
      os << *itForbClassTrain << " " << std::endl;
    }   
    os << this->createEndTag( "forbidden_classesTrain" ) << std::endl;
    
    //
    
    os << this->createStartTag( "forbidden_classesActiveLearning" ) << std::endl;
    os << "size: " << forbidden_classesActiveLearning.size() << std::endl;

    for ( std::set< int >::const_iterator itForbClassAL = forbidden_classesActiveLearning.begin();
          itForbClassAL != forbidden_classesActiveLearning.end();
          itForbClassAL++
        )
    {
      os << *itForbClassAL << " " << std::endl;
    }   
    os << this->createEndTag( "forbidden_classesActiveLearning" ) << std::endl;     
    
    //
    
    os << this->createStartTag( "classesInUse" ) << std::endl;
    os << "size: " << classesInUse.size() << std::endl;

    for ( std::set< int >::const_iterator itClassesInUse = classesInUse.begin();
          itClassesInUse != classesInUse.end();
          itClassesInUse++
        )
    {
      os << *itClassesInUse << " " << std::endl;
    }   
    os << this->createEndTag( "classesInUse" ) << std::endl;
    
    
    os << this->createStartTag( "numberOfClasses" ) << std::endl;
    os << this->numberOfClasses << std::endl;
    os << this->createStartTag( "numberOfClasses" ) << std::endl;     
    
    
    os << this->createStartTag( "read_classifier" ) << std::endl;
    os << this->read_classifier << std::endl;
    os << this->createStartTag( "read_classifier" ) << std::endl;     
    
    
    os << this->createStartTag( "save_classifier" ) << std::endl;
    os << this->save_classifier << std::endl;
    os << this->createStartTag( "save_classifier" ) << std::endl;     
    
    
    os << this->createStartTag( "cache" ) << std::endl;
    os << this->cache << std::endl;
    os << this->createStartTag( "cache" ) << std::endl;     
    
    
    os << this->createStartTag( "resultdir" ) << std::endl;
    os << this->resultdir << std::endl;
    os << this->createStartTag( "resultdir" ) << std::endl;     

    //TODO newTrainExamples    
    
    ///////////////////////////////
    //     SEGMENTATION STUFF    //
    ///////////////////////////////
    
    // theoretically, we should properly store and restore the regionSeg object. However,  its parent class does not provide
    // a Persistent interface yet. Hence, we perform this tiny workaround which works, since regionSeg is not changed over time...
    // only be aware of parameters originally set via config...    
    os << this->createStartTag( "s_rsMethode" ) << std::endl;
    os << this->s_rsMethode << std::endl;
    os << this->createStartTag( "s_rsMethode" ) << std::endl;    
    
    os << this->createStartTag( "reuseSegmentation" ) << std::endl;
    os << this->reuseSegmentation << std::endl;
    os << this->createStartTag( "reuseSegmentation" ) << std::endl;    
    
    os << this->createStartTag( "queriedRegions" ) << std::endl;
    os << "size: " << queriedRegions.size() << std::endl;
    std::map< std::string, std::set< int > >::const_iterator itQueriedRegions = queriedRegions.begin();
    for ( uint i = 0; i < queriedRegions.size(); i++ )
    {
      // store key
      os << itQueriedRegions->first << std::endl;
      
      // store values -- inner loop over sets
      os << "size: " << ( itQueriedRegions->second ).size() << std::endl;

      for ( std::set< int >::const_iterator itRegionsOfImg = ( itQueriedRegions->second ).begin();
	    itRegionsOfImg != ( itQueriedRegions->second ).end();
	    itRegionsOfImg++
	  )
      {
	os << *itRegionsOfImg << " " << std::endl;
      }       
      
      itQueriedRegions++;
    } 
    os << this->createStartTag( "queriedRegions" ) << std::endl;
    //
    //TODO currentRegionToQuery

    
    
    ///////////////////////////////
    //       PARENT OBJECT       //
    ///////////////////////////////
    os << this->createStartTag( "SemSegNovelty--Parent" ) << std::endl;
    SemanticSegmentation::store(os);  
    os << this->createStartTag( "SemSegNovelty--Parent" ) << std::endl;
    
    
    // done
    os << this->createEndTag( "SemSegNovelty" ) << std::endl;    
  }
  else
  {
    std::cerr << "OutStream not initialized - storing not possible!" << std::endl;
  }
}

void SemSegNovelty::clear ()
{
 //TODO
}