#include <sstream>
#include <iostream>

#include "SemSegCsurka.h"
#include "vislearning/baselib/ICETools.h"
#include "core/image/FilterT.h"
#include "semseg/semseg/postsegmentation/PSSImageLevelPrior.h"

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

#undef DEBUG_CSURK

#undef UNCERTAINTY
// #define UNCERTAINTY

SemSegCsurka::SemSegCsurka ( const Config *conf,
                             const MultiDataset *md )
    : SemanticSegmentation ( conf, & ( md->getClassNames ( "train" ) ) )
{
  this->conf = conf;

  opSiftImpl = conf->gS ( "Descriptor", "implementation", "VANDESANDE" );
  readfeat = conf->gB ( "Descriptor", "read", true );
  writefeat = conf->gB ( "Descriptor", "write", true );
#ifdef DEBUG_CSURK
  clog << "[log] SemSegCsurka::SemSegCsurka: OppenentSift implemenation: " << opSiftImpl << endl;
#endif

  save_cache = conf->gB ( "FPCPixel", "save_cache", true );
  read_cache = conf->gB ( "FPCPixel", "read_cache", false );
  cache = conf->gS ( "cache", "root", "" );
  sigmaweight = conf->gD ( "SemSegCsurka", "sigmaweight", 0.6 );

  dim = conf->gI ( "SemSegCsurka", "pcadim", 50 );

  usepca = conf->gB ( "SemSegCsurka", "usepca", true );
  calcpca = conf->gB ( "SemSegCsurka", "calcpca", false );

  usegmm = conf->gB ( "SemSegCsurka", "usegmm", false );
  norm = conf->gB ( "SemSegCsurka", "normalize", false );
  usefisher = conf->gB ( "SemSegCsurka", "usefisher", false );
  dogmm = conf->gB ( "SemSegCsurka", "dogmm", false );
  gaussians = conf->gI ( "SemSegCsurka", "gaussians", 50 );

  usekmeans = conf->gB ( "SemSegCsurka", "usekmeans", false );
  kmeansfeat = conf->gI ( "SemSegCsurka", "kmeansfeat", 50 );
  kmeanshard = conf->gB ( "SemSegCsurka", "kmeanshard", false );

  cname = conf->gS ( "SemSegCsurka", "classifier", "RandomForests" );
  anteil = conf->gD ( "SemSegCsurka", "anteil", 1.0 );
  userellocprior = conf->gB ( "SemSegCsurka", "rellocfeat", false );
  bool usesrg = conf->gB ( "SemSegCsurka", "usesrg", false );

  useregions = conf->gB ( "SemSegCsurka", "useregions", true );
  savesteps = conf->gB ( "SemSegCsurka", "savesteps", true );
  bool usegcopt = conf->gB ( "SemSegCsurka", "usegcopt", false );

  bestclasses = conf->gI ( "SemSegCsurka", "bestclasses", 0 );

  smoothhl = conf->gB ( "SemSegCsurka", "smoothhl", false );
  smoothfactor = conf->gD ( "SemSegCsurka", "smoothfactor", 1.0 );

  usecolorfeats = conf->gB ( "SemSegCsurka", "usecolorfeats", false );

  string rsMethod = conf->gS ( "SemSegCsurka", "segmentation", "meanshift" );

  g = NULL;
  k = NULL;
  relloc = NULL;
  srg = NULL;
  gcopt = NULL;

  if ( !useregions && ( userellocprior || usesrg ) )
  {
    cerr << "relative location priors and super region growing are just supported in combination with useregions" << endl;
    exit ( 1 );
  }

  if ( usepca )
    pca = PCA ( dim );

  RegionSegmentationMethod * tmpseg;
  if ( rsMethod == "meanshift" )
    tmpseg = new RSMeanShift ( conf );
  else
    tmpseg = new RSGraphBased ( conf );

  if ( save_cache )
    seg = new RSCache ( conf, tmpseg );
  else
    seg = tmpseg;

  if ( userellocprior )
    relloc = new RelativeLocationPrior ( conf );
  else
    relloc = NULL;

#ifdef NICE_USELIB_ICE
  if ( usesrg )
    srg = new PPSuperregion ( conf );
  else
    srg = NULL;
#else
  srg = NULL;
#endif

  if ( usegcopt )
    gcopt = new PPGraphCut ( conf );
  else
    gcopt = NULL;

  classifier = NULL;
  vclassifier = NULL;
  if ( cname == "RandomForests" )
    classifier = new FPCRandomForests ( conf, "ClassifierForest" );
  else if ( cname == "SMLR" )
    classifier = new FPCSMLR ( conf, "ClassifierSMLR" );
  else if ( cname == "GPHIK" )
    classifier = new FPCGPHIK ( conf, "ClassiferGPHIK" );
  else
    vclassifier = GenericClassifierSelection::selectVecClassifier ( conf, "main" );
  //classifier = new FPCSparseMultinomialLogisticRegression(conf, "ClassifierSMLR");

  if ( classifier != NULL )
    classifier->setMaxClassNo ( classNames->getMaxClassno() );
  else
    vclassifier->setMaxClassNo ( classNames->getMaxClassno() );

  cn = md->getClassNames ( "train" );

  if ( read_cache )
  {
    fprintf ( stderr, "SemSegCsurka:: Reading classifier data from %s\n", ( cache + "/fpcrf.data" ).c_str() );

    if ( classifier != NULL )
      classifier->read ( cache + "/fpcrf.data" );
    else
      vclassifier->read ( cache + "/veccl.data" );

    if ( usepca )
    {
      std::string filename = cache + "/pca";
      pca.read ( filename );
    }

    if ( usegmm )
    {
      g = new GMM ( conf, gaussians );

      if ( !g->loadData ( cache + "/gmm" ) )
      {
        cerr << "SemSegCsurka:: no gmm file found" << endl;
        exit ( -1 );
      }
    }
    else {
      g = NULL;
    }

    if ( usekmeans )
    {
      k = new KMeansOnline ( gaussians );
    }

    fprintf ( stderr, "SemSegCsurka:: successfully read\n" );

    std::string filename = cache + "/rlp";

    FILE *value;
    value = fopen ( filename.c_str(), "r" );

    if ( value == NULL )
    {
      trainpostprocess ( md );
    }
    else
    {
      if ( userellocprior )
      {
        relloc->read ( filename );
      }
    }

    filename = cache + "/srg";

    value = fopen ( filename.c_str(), "r" );

    if ( value == NULL )
    {
      trainpostprocess ( md );
    }
    else
    {
      if ( srg != NULL )
      {
        srg->read ( filename );
      }
    }
  }
  else
  {
    train ( md );
  }
}

SemSegCsurka::~SemSegCsurka()
{
  // clean-up
  if ( classifier != NULL )
    delete classifier;
  if ( vclassifier != NULL )
    delete vclassifier;
  if ( seg != NULL )
    delete seg;

  g = NULL;
  if ( g != NULL )
    delete g;
}

void SemSegCsurka::normalize ( Examples &ex )
{
  assert ( ex.size() > 0 );
  if ( vecmin.size() == 0 )
  {
    for ( int j = 0; j < ( int ) ex[0].second.vec->size(); j++ )
    {
      double maxv = -numeric_limits<int>::max();
      double minv = numeric_limits<int>::max();
      for ( int i = 0; i < ( int ) ex.size(); i++ )
      {
        maxv = std::max ( maxv, ( *ex[i].second.vec ) [j] );
        minv = std::min ( minv, ( *ex[i].second.vec ) [j] );
      }
      vecmin.push_back ( minv );
      vecmax.push_back ( maxv );
    }
  }
  for ( int i = 0; i < ( int ) ex.size(); i++ )
  {
    for ( int j = 0; j < ( int ) ex[i].second.vec->size(); j++ )
    {
      ( *ex[i].second.vec ) [j] = ( ( *ex[i].second.vec ) [j] - vecmin[j] ) / ( vecmax[j] - vecmin[j] );
    }
  }
  return;
}

void SemSegCsurka::convertLowToHigh ( Examples &ex, double reduce )
{
  cout << "converting low-level features to high-level features" << endl;

  if ( reduce >= 1.0 )
  {
    for ( int i = 0; i < ( int ) ex.size(); i++ )
    {
      SparseVector *f = new SparseVector();

      if ( usekmeans )
      {
        k->getDist ( *ex[i].second.vec, *f, kmeansfeat, kmeanshard );
      }
      else
      {
        if ( usefisher )
          g->getFisher ( *ex[i].second.vec, *f );
        else
          g->getProbs ( *ex[i].second.vec, *f );
      }
      delete ex[i].second.vec;

      ex[i].second.vec = NULL;
      ex[i].second.svec = f;
    }
  }
  else
  {
    srand ( time ( NULL ) );

    vector<bool> del ( ex.size(), false );
    cout << "Example size old " << ex.size() << endl;

//#pragma omp parallel for
    for ( int i = 0; i < ( int ) ex.size(); i++ )
    {
      double rval = ( double ) rand() / ( double ) RAND_MAX;
      if ( rval < reduce )
      {
        SparseVector *f = new SparseVector();

        if ( usekmeans )
          k->getDist ( *ex[i].second.vec, *f, kmeansfeat, kmeanshard );
        else
        {
          if ( usefisher )
            g->getFisher ( *ex[i].second.vec, *f );
          else
            g->getProbs ( *ex[i].second.vec, *f );
        }

        delete ex[i].second.vec;
        ex[i].second.vec = NULL;
        ex[i].second.svec = f;
      }
      else
      {
        del[i] = true;
      }
    }
    for ( int i = ( int ) del.size() - 1; i >= 0; i-- )
    {
      if ( del[i] )
      {
        ex.erase ( ex.begin() + i );
      }
    }
    cerr << "Example size new " << ex.size() << endl;
  }
  cerr << "converting low-level features to high-level features finished" << endl;
}

void SemSegCsurka::smoothHL ( Examples ex )
{

  if ( !smoothhl )
    return;
  assert ( ex.size() > 1 );

  long long int minx = numeric_limits<long long int>::max();
  long long int miny = numeric_limits<long long int>::max();
  long long int maxx = -numeric_limits<long long int>::max();
  long long int maxy = -numeric_limits<long long int>::max();
  long long int distx = numeric_limits<long long int>::max();
  long long int disty = numeric_limits<long long int>::max();

  set<double> scales;
  for ( int i = 0; i < ( int ) ex.size(); i++ )
  {
    scales.insert ( ex[i].second.scale );
  }

  map<double, int> scalepos;
  int it = 0;

  for ( set<double>::const_iterator iter = scales.begin(); iter != scales.end();    ++iter, ++it )
  {
    scalepos.insert ( make_pair ( *iter, it ) );
  }

  for ( int i = 0; i < ( int ) ex.size(); i++ )
  {
    if ( minx < numeric_limits<int>::max() && ex[i].second.x - minx > 0 )
      distx = std::min ( distx, ex[i].second.x - minx );
    if ( miny < numeric_limits<int>::max() && ex[i].second.y - miny > 0 )
      disty = std::min ( disty, ex[i].second.y - miny );
    minx = std::min ( ( long long int ) ex[i].second.x, minx );
    maxx = std::max ( ( long long int ) ex[i].second.x, maxx );
    miny = std::min ( ( long long int ) ex[i].second.y, miny );
    maxy = std::max ( ( long long int ) ex[i].second.y, maxy );
  }

  distx = abs ( distx );

  int xsize = ( maxx - minx ) / distx + 1;
  int ysize = ( maxy - miny ) / disty + 1;
  double valx = ( ( double ) xsize - 1 ) / ( double ) ( maxx - minx );
  double valy = ( ( double ) ysize - 1 ) / ( double ) ( maxy - miny );

  //double sigma = smoothfactor;
  double sigma = std::max ( xsize, ysize ) * smoothfactor;
  //double sigma = 0.2;
  cout << "sigma1: " << sigma << endl;

  vector<NICE::FloatImage> imgv;
  vector<NICE::FloatImage> gaussImgv;
  for ( int i = 0; i < ( int ) scalepos.size(); i++ )
  {
    NICE::FloatImage img ( xsize, ysize );
    NICE::FloatImage gaussImg ( xsize, ysize );
    imgv.push_back ( img );
    gaussImgv.push_back ( gaussImg );
  }

  for ( int d = 0; d < ex[0].second.svec->getDim(); d++ )
  {
    //TODO: max und min dynamisches bestimmen

    for ( int i = 0; i < ( int ) scalepos.size(); i++ )
    {
      imgv[i].set ( 0.0 );
      gaussImgv[i].set ( 0.0 );
    }

    for ( int i = 0; i < ( int ) ex.size(); i++ )
    {
      int xpos = ( ex[i].second.x - minx ) * valx;
      int ypos = ( ex[i].second.y - miny ) * valy;

      double val = ex[i].second.svec->get ( d );
      imgv[scalepos[ex[i].second.scale]].setPixel ( xpos, ypos, val );
    }

    /*
    for(int y = 0; y < ysize; y++)
    {
     for(int x = 0; x < xsize; x++)
     {
      // refactor-nice.pl: check this substitution
      // old: double val = GetValD(img,x,y);
      double val = img.getPixel(x,y);
      double  c = 0.0;
      if(val == 0.0)
      {
       if(x > 0)
       {
        // refactor-nice.pl: check this substitution
        // old: val+=GetValD(img,x-1,y);
        val+=img.getPixel(x-1,y);
        c+=1.0;
       }
       if(y > 0)
       {
        // refactor-nice.pl: check this substitution
        // old: val+=GetValD(img,x,y-1);
        val+=img.getPixel(x,y-1);
        c+=1.0;
       }
       if(x < xsize-1)
       {
        // refactor-nice.pl: check this substitution
        // old: val+=GetValD(img,x+1,y);
        val+=img.getPixel(x+1,y);
        c+=1.0;
       }
       if(y < ysize-1)
       {
        // refactor-nice.pl: check this substitution
        // old: val+=GetValD(img,x,y+1);
        val+=img.getPixel(x,y+1);
        c+=1.0;
       }
       // refactor-nice.pl: check this substitution
       // old: PutValD(img,x,y,val/c);
       img.setPixel(x,y,val/c);
      }
     }
    }*/

    for ( int i = 0; i < ( int ) imgv.size(); i++ )
    {
      NICE::FilterT<float, float, float> filt;
      filt.filterGaussSigmaApproximate ( imgv[i], sigma, &gaussImgv[i] );
    }

    for ( int i = 0; i < ( int ) ex.size(); i++ )
    {
      int xpos = ( ex[i].second.x - minx ) * valx;
      int ypos = ( ex[i].second.y - miny ) * valy;
      // refactor-nice.pl: check this substitution
      // old: double val = GetValD ( gaussImgv[scalepos[ex[i].second.scale]], xpos, ypos );
      double val = gaussImgv[scalepos[ex[i].second.scale]].getPixel ( xpos, ypos );

      if ( fabs ( val ) < 1e-7 )
      {
        if ( ex[i].second.svec->get ( d ) != 0.0 )
        {
          ex[i].second.svec->erase ( d );
        }
      }
      else
      {
        ( *ex[i].second.svec ) [d] = val;
      }
    }
  }
}

void SemSegCsurka::initializePCA ( Examples &ex )
{
#ifdef DEBUG
  cerr << "start computing pca" << endl;
#endif
  std::string filename = cache + "/pca";
  FILE *value;
  value = fopen ( filename.c_str(), "r" );

  if ( value == NULL || calcpca )
  {
    srand ( time ( NULL ) );

    int featsize = ( int ) ex.size();
    int maxfeatures = dim * 10;
    int olddim = ex[0].second.vec->size();

    maxfeatures = std::min ( maxfeatures, featsize );

    NICE::Matrix features ( maxfeatures, olddim );

    for ( int i = 0; i < maxfeatures; i++ )
    {
      int k = rand() % featsize;

      int vsize = ( int ) ex[k].second.vec->size();
      for ( int j = 0; j < vsize; j++ )
      {
        features ( i, j ) = ( * ( ex[k].second.vec ) ) [j];
      }
    }
    pca.calculateBasis ( features, dim );

    if ( save_cache )
      pca.save ( filename );

  }
  else
  {
    cout << "readpca: " << filename << endl;
    pca.read ( filename );
    cout << "end" << endl;
  }
#ifdef DEBUG
  cerr << "finished computing pca" << endl;
#endif
}

void SemSegCsurka::doPCA ( Examples &ex )
{
  cout << "converting features using pca starts" << endl;

  std::string savedir = cname = conf->gS ( "cache", "root", "/dev/null/" );
  std::string shortf = ex.filename;
  if ( string::npos != ex.filename.rfind ( "/" ) )
    shortf = ex.filename.substr ( ex.filename.rfind ( "/" ) );
  std::string filename = savedir + "/pcasave/" + shortf;
  std::string syscall = "mkdir " + savedir + "/pcasave";
  system ( syscall.c_str() );
  cout << "filename: " << filename << endl;

  if ( !FileMgt::fileExists ( filename ) || calcpca )
  {
    ofstream ofStream;

    //Opens the file binary
    ofStream.open ( filename.c_str(), fstream::out | fstream::binary );

    for ( int k = 0; k < ( int ) ex.size(); k++ )
    {
      NICE::Vector tmp = pca.getFeatureVector ( * ( ex[k].second.vec ), true );
      delete ex[k].second.vec;
      for ( int d = 0; d < ( int ) tmp.size(); d++ )
        ofStream.write ( ( char* ) &tmp[d], sizeof ( double ) );
      ex[k].second.vec = new NICE::Vector ( tmp );
    }
    ofStream.close();
    cout << endl;
  }
  else
  {
    ifstream ifStream;
    ifStream.open ( filename.c_str(), std::fstream::in | std::fstream::binary );
    for ( int k = 0; k < ( int ) ex.size(); k++ )
    {
      NICE::Vector tmp = NICE::Vector ( dim );
      delete ex[k].second.vec;
      for ( int d = 0; d < dim; d++ )
        ifStream.read ( ( char* ) &tmp[d], sizeof ( double ) );
      ex[k].second.vec = new NICE::Vector ( tmp );
    }

    ifStream.close();
  }
  cout << "converting features using pca finished" << endl;
}

void SemSegCsurka::train ( const MultiDataset *md )
{

  /*die einzelnen Trainingsschritte
  1. auf allen Trainingsbilder SIFT Merkmale an den Gitterpunkten bei allen Auflösungen bestimmen
  2. PCA anwenden
  3. aus diesen ein GMM erstellen
  4. für jedes SIFT-Merkmal einen Vektor erstellen, der an der Stelle i die Wahrscheinlichkeit enthällt zur Verteilung i des GMM, Zur Zeit mit BoV-Alternative durch Moosman06 erledigt
  5. diese Vektoren in einem diskriminitativen Klassifikator ( z.B. SLR oder Randomized Forests) zusammen mit ihrer Klassenzugehörigkeit anlernen
  */
#ifdef DEBUG
  cerr << "SemSegCsurka:: training starts" << endl;
#endif

  Examples examples;
  examples.filename = "training";


  // Welche Opponentsift Implementierung soll genutzt werden ?
  LocalFeatureRepresentation *cSIFT = NULL;
  LocalFeatureRepresentation *writeFeats = NULL;
  LocalFeatureRepresentation *readFeats = NULL;
  LocalFeatureRepresentation *getFeats = NULL;

  if ( opSiftImpl == "NICE" )
  {
    cSIFT = new LFonHSG ( conf, "HSGtrain" );
  }
  else if ( opSiftImpl == "VANDESANDE" )
  {
    // the used features
    cSIFT = new LFColorSande ( conf, "LFColorSandeTrain" );
  }
  else
  {
    fthrow ( Exception, "feattype: %s not yet supported" << opSiftImpl );
  }

  getFeats = cSIFT;

  if ( writefeat )
  {
    // write the features to a file, if there isn't any to read
    writeFeats = new LFWriteCache ( conf, cSIFT );
    getFeats = writeFeats;
  }

  if ( readfeat )
  {
    // read the features from a file
    if ( writefeat )
    {
      readFeats = new LFReadCache ( conf, writeFeats, -1 );
    }
    else
    {
      readFeats = new LFReadCache ( conf, cSIFT, -1 );
    }
    getFeats = readFeats;
  }

  // additional Colorfeatures
  LocalFeatureColorWeijer lcw ( conf );

  int lfdimension = -1;

  const LabeledSet train = * ( *md ) ["train"];
  const LabeledSet *trainp = &train;

  ////////////////////////
  // Merkmale berechnen //
  ////////////////////////

  std::string forbidden_classes_s = conf->gS ( "analysis", "donttrain", "" );
  if ( forbidden_classes_s == "" )
  {
    forbidden_classes_s = conf->gS ( "analysis", "forbidden_classes", "" );
  }
  cn.getSelection ( forbidden_classes_s, forbidden_classes );
  cerr << "forbidden: " << forbidden_classes_s << endl;

  ProgressBar pb ( "Local Feature Extraction" );
  pb.show();

  int imgnb = 0;

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

    pb.update ( trainp->count() );

    NICE::ColorImage img;

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

    fprintf ( stderr, "SemSegCsurka: Collecting pixel examples from localization info: %s\n",
              currentFile.c_str() );

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

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

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

    Globals::setCurrentImgFN ( currentFile );

    VVector features;
    VVector cfeatures;
    VVector positions;

    NICE::ColorImage cimg ( currentFile );

    getFeats->extractFeatures ( img, features, positions );

#ifdef DEBUG_CSURK
    cout << "[log] SemSegCsruka::train -> " << currentFile << " an " << positions.size() << " Positionen wurden Features (Anz = " << features.size() << ") " << endl;
    cout << "mit einer Dimension von " << features[ 0].size() << " extrahiert." << endl;
#endif

    if ( usecolorfeats )
      lcw.getDescriptors ( cimg, cfeatures, positions );

    int j = 0;

    for ( VVector::const_iterator i = features.begin();
          i != features.end();
          i++, j++ )
    {
      const NICE::Vector & x = *i;
      classno = pixelLabels.getPixel ( ( int ) positions[j][0], ( int ) positions[j][1] );

      if ( forbidden_classes.find ( classno ) != forbidden_classes.end() )
        continue;

      if ( lfdimension < 0 )
        lfdimension = ( int ) x.size();
      else
        assert ( lfdimension == ( int ) x.size() );

      NICE::Vector *v = new NICE::Vector ( x );

      if ( usecolorfeats && !usepca )
        v->append ( cfeatures[j] );

      Example example ( v );
      example.position = imgnb;
      examples.push_back (
        pair<int, Example> ( classno, example ) );
    }
    features.clear();
    positions.clear();
    delete ce;
    imgnb++;
  }

  pb.hide();

  //////////////////
  // PCA anwenden //
  //////////////////

  if ( usepca )
  {
    if ( !read_cache )
    {
      initializePCA ( examples );
    }
    doPCA ( examples );
    lfdimension = dim;
  }

  /////////////////////////////////////////////////////
  // Low-Level Features in High-Level transformieren //
  /////////////////////////////////////////////////////

  int hlfdimension = lfdimension;

  if ( norm )
    normalize ( examples );

  if ( usegmm )
  {
    if ( !usepca && !norm )
      normalize ( examples );
    g = new GMM ( conf, gaussians );

    if ( dogmm || !g->loadData ( cache + "/gmm" ) )
    {
      g->computeMixture ( examples );
      if ( save_cache )
        g->saveData ( cache + "/gmm" );
    }

    hlfdimension = gaussians;

    if ( usefisher )
      hlfdimension = gaussians * 2 * dim;
  }

  if ( usekmeans )
  {
    if ( !usepca || norm )
      normalize ( examples );
    k = new KMeansOnline ( gaussians );

    k->cluster ( examples );

    hlfdimension = gaussians;
  }

  if ( usekmeans || usegmm )
  {
    examples.clear();
    pb.reset ( "Local Feature Extraction" );
    lfdimension = -1;
    pb.update ( trainp->count() );
    LOOP_ALL_S ( *trainp )
    {
      EACH_INFO ( classno, info );

      pb.update ( trainp->count() );

      NICE::ColorImage img;

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

      fprintf ( stderr, "SemSegCsurka: Collecting pixel examples from localization info: %s\n",
                currentFile.c_str() );

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

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

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

      Globals::setCurrentImgFN ( currentFile );

      VVector features;
      VVector cfeatures;
      VVector positions;

      NICE::ColorImage cimg ( currentFile );

      getFeats->extractFeatures ( img, features, positions );

      if ( usecolorfeats )
        lcw.getDescriptors ( cimg, cfeatures, positions );

      int j = 0;

      Examples tmpex;

      for ( VVector::const_iterator i = features.begin();
            i != features.end();
            i++, j++ )
      {

        const NICE::Vector & x = *i;

        classno = pixelLabels.getPixel ( ( int ) positions[j][0], ( int ) positions[j][1] );

        if ( forbidden_classes.find ( classno ) != forbidden_classes.end() )
          continue;

        if ( lfdimension < 0 )
          lfdimension = ( int ) x.size();
        else
          assert ( lfdimension == ( int ) x.size() );

        NICE::Vector *v = new NICE::Vector ( x );
        if ( usecolorfeats )
          v->append ( cfeatures[j] );

        Example example ( v );
        example.position = imgnb;
        example.x = ( int ) positions[j][0];
        example.y = ( int ) positions[j][1];
        example.scale = positions[j][2];

        tmpex.push_back ( pair<int, Example> ( classno, example ) );
      }
      tmpex.filename = currentFile;
      if ( usepca )
      {
        doPCA ( tmpex );
      }

      convertLowToHigh ( tmpex, anteil );

      smoothHL ( tmpex );

      for ( int i = 0; i < ( int ) tmpex.size(); i++ )
      {
        examples.push_back ( pair<int, Example> ( tmpex[i].first, tmpex[i].second ) );
      }

      tmpex.clear();

      features.clear();
      positions.clear();
      delete ce;
      imgnb++;

    }

    pb.hide();
  }
  ////////////////////////////
  // Klassifikator anlernen //
  ////////////////////////////
  FeaturePool fp;

  Feature *f;

  if ( usegmm || usekmeans )
    f = new SparseVectorFeature ( hlfdimension );
  else
    f = new VectorFeature ( hlfdimension );

  f->explode ( fp );
  delete f;

  if ( usecolorfeats && ! ( usekmeans || usegmm ) )
  {
    int dimension = hlfdimension + 11;
    for ( int i = hlfdimension ; i < dimension ; i++ )
    {
      VectorFeature *f = new VectorFeature ( dimension );
      f->feature_index = i;
      fp.addFeature ( f, 1.0 / dimension );
    }
  }
  /*
  cout << "train classifier" << endl;
  fp.store(cout);
  getchar();
  for(int z = 0; z < examples.size(); z++)
  {
  cout << "examples.size() " << examples.size() << endl;
  cout << "class: " << examples[z].first << endl;
   cout << *examples[z].second.vec << endl;
   getchar();
  }*/

  if ( classifier != NULL )
    classifier->train ( fp, examples );
  else
  {
    LabeledSetVector lvec;
    convertExamplesToLSet ( examples, lvec );
    vclassifier->teach ( lvec );
    if ( usegmm )
      convertLSetToSparseExamples ( examples, lvec );
    else
      convertLSetToExamples ( examples, lvec );
    vclassifier->finishTeaching();
  }

  fp.destroy();

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

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

  if ( cSIFT != NULL )
    delete cSIFT;
  if ( writeFeats != NULL )
    delete writeFeats;
  if ( readFeats != NULL )
    delete readFeats;
  getFeats = NULL;

  trainpostprocess ( md );

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

void SemSegCsurka::trainpostprocess ( const MultiDataset *md )
{
  cout << "start postprocess" << endl;
  ////////////////////////////
  // Postprocess trainieren //
  ////////////////////////////
  const LabeledSet train = * ( *md ) ["train"];
  const LabeledSet *trainp = &train;

  if ( userellocprior || srg != NULL || gcopt != NULL )
  {
    clog << "[log] SemSegCsurka::trainpostprocess: if ( userellocprior || srg != NULL || gcopt !=NULL )" << endl;
    if ( userellocprior )
      relloc->setClassNo ( cn.numClasses() );

    if ( gcopt != NULL )
    {
      gcopt->setClassNo ( cn.numClasses() );
    }

    ProgressBar pb ( "learn relative location prior maps" );
    pb.show();
    LOOP_ALL_S ( *trainp ) // für alle Bilder den ersten Klassifikationsschritt durchführen um den zweiten Klassifikator anzutrainieren
    {
      EACH_INFO ( classno, info );

      pb.update ( trainp->count() );

      NICE::ColorImage img;

      std::string currentFile = info.img();
      Globals::setCurrentImgFN ( currentFile );
      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;
      }

      fprintf ( stderr, "SemSegCsurka: Collecting pixel examples from localization info: %s\n",
                currentFile.c_str() );

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

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

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

      //Regionen ermitteln
      NICE::Matrix mask;

      int regionsize = seg->segRegions ( img, mask );
#ifdef DEBUG_CSURK
      Image overlay ( img.width(), img.height() );

      double maxval = -numeric_limits<double>::max();

      for ( int y = 0; y < img.height(); y++ )
      {
        for ( int x = 0; x < img.width(); x++ )
        {
          int val = ( ( int ) mask ( x, y ) + 1 ) % 256;
          overlay.setPixel ( x, y, val );
          maxval = std::max ( mask ( x, y ), maxval );
        }
      }

      cout << maxval << " different regions found" << endl;

      NICE::showImageOverlay ( img, overlay, "Segmentation Result" );
#endif

      Examples regions;

      vector<vector<int> > hists;

      for ( int i = 0; i < regionsize; i++ )
      {
        Example tmp;
        regions.push_back ( pair<int, Example> ( 0, tmp ) );
        vector<int> hist ( cn.numClasses(), 0 );
        hists.push_back ( hist );
      }

      for ( int x = 0; x < xsize; x++ )
      {
        for ( int y = 0; y < ysize; y++ )
        {
          int numb = mask ( x, y );
          regions[numb].second.x += x;
          regions[numb].second.y += y;
          regions[numb].second.weight += 1.0;
          hists[numb][pixelLabels.getPixel ( x,y ) ]++;
        }
      }

      for ( int i = 0; i < regionsize; i++ )
      {
        regions[i].second.x /= ( int ) regions[i].second.weight;
        regions[i].second.y /= ( int ) regions[i].second.weight;

        int maxval = -numeric_limits<int>::max();
        int maxpos = -1;
        int secondpos = -1;
        for ( int k = 0; k < ( int ) hists[i].size(); k++ )
        {
          if ( maxval < hists[i][k] )
          {
            maxval = hists[i][k];
            secondpos = maxpos;
            maxpos = k;
          }
        }

        if ( cn.text ( maxpos ) == "various" )
          regions[i].first = secondpos;
        else
          regions[i].first = maxpos;

      }
      if ( userellocprior )
        relloc->trainPriorsMaps ( regions, xsize, ysize );

      if ( srg != NULL )
        srg->trainShape ( regions, mask );

      if ( gcopt != NULL )
        gcopt->trainImage ( regions, mask );

      delete ce;

    }
    pb.hide();
    if ( userellocprior )
      relloc->finishPriorsMaps ( cn );

    if ( srg != NULL )
      srg->finishShape ( cn );

    if ( gcopt != NULL )
      gcopt->finishPP ( cn );
  }
  if ( userellocprior )
  {
    clog << "[log] SemSegCsurka::trainpostprocess: if ( userellocprior )" << endl;
    ProgressBar pb ( "learn relative location classifier" );
    pb.show();

    int nummer = 0;
    LOOP_ALL_S ( *trainp ) // für alle Bilder den ersten Klassifikationsschritt durchführen um den zweiten Klassifikator anzutrainieren
    {
      //EACH_S(classno, currentFile);
      EACH_INFO ( classno, info );
      nummer++;
      pb.update ( trainp->count() );

      NICE::Image img;
      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;
      }

      fprintf ( stderr, "SemSegCsurka: Collecting pixel examples from localization info: %s\n",
                currentFile.c_str() );

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

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

      try {
        img = Preprocess::ReadImgAdv ( currentFile.c_str() );
      }
      catch ( Exception )
      {
        cerr << "SemSegCsurka: error opening image file <" << currentFile << ">" << endl;
        continue;
      }
      Globals::setCurrentImgFN ( currentFile );

      NICE::ImageT<int> segresult;

      NICE::MultiChannelImageT<double> probabilities ( xsize, ysize, classno );

      Examples regions;

      NICE::Matrix mask;

      if ( savesteps )
      {
        std::ostringstream s1;
        s1 << cache << "/rlpsave/" << nummer;

        std::string filename = s1.str();
        s1 << ".probs";

        std::string fn2 = s1.str();

        FILE *file;
        file = fopen ( filename.c_str(), "r" );

        if ( file == NULL )
        {
          //berechnen
          classifyregions ( ce, segresult, probabilities, regions, mask );
          //schreiben
          ofstream fout ( filename.c_str(), ios::app );
          fout << regions.size() << endl;
          for ( int i = 0; i < ( int ) regions.size(); i++ )
          {
            regions[i].second.store ( fout );
            fout << regions[i].first << endl;
          }
          fout.close();
          probabilities.store ( fn2 );
        }
        else
        {
          //lesen
          ifstream fin ( filename.c_str() );
          int size;
          fin >> size;

          for ( int i = 0; i < size; i++ )
          {
            Example ex;
            ex.restore ( fin );
            int tmp;
            fin >> tmp;
            regions.push_back ( pair<int, Example> ( tmp, ex ) );
          }

          fin.close();

          probabilities.restore ( fn2 );
        }
      }
      else
      {
        classifyregions ( ce, segresult, probabilities, regions, mask );
      }

      relloc->trainClassifier ( regions, probabilities );

      delete ce;

    }
    relloc->finishClassifier();
    pb.hide();

    relloc->save ( cache + "/rlp" );
  }
  cout << "finished postprocess" << endl;
}

void SemSegCsurka::classifyregions ( CachedExample *ce, NICE::ImageT<int> & segresult, NICE::MultiChannelImageT<double> & probabilities, Examples &Regionen, NICE::Matrix & mask )
{
  /* die einzelnen Testschritte:
  1.x  auf dem Testbild alle SIFT Merkmale an den Gitterpunkten bei allen Auflösungen bestimmen
  2.x  für jedes SIFT-Merkmal einen Vektor erstellen, der an der Stelle i die Wahrscheinlichkeit enthällt zur Verteilung i des GMM
  3.x diese Vektoren klassifizieren, so dass für jede Klasse die Wahrscheinlichkeit gespeichert wird
  4.x für jeden Pixel die Wahrscheinlichkeiten mitteln aus allen Patches, in denen der Pixel vorkommt
  5.x das Originalbild in homogene Bereiche segmentieren
  6.x die homogenen Bereiche bekommen die gemittelten Wahrscheinlichkeiten ihrer Pixel
  7. (einzelne Klassen mit einem globalen Klassifikator ausschließen)
  8.x jeder Pixel bekommt die Klasse seiner Region zugeordnet
  */

  clog << "[log] SemSegCsruka::classifyregions" << endl;
  int xsize, ysize;

  ce->getImageSize ( xsize, ysize );

  probabilities.reInit ( xsize, ysize, classNames->getMaxClassno() + 1 );
  clog << "[log] SemSegCsruka::classifyregions: probabilities.channels() = " << probabilities.channels() << endl;

  segresult.resize ( xsize, ysize );

  Examples pce;

  // Welche Opponentsift Implementierung soll genutzt werden ?
  LocalFeatureRepresentation *cSIFT = NULL;
  LocalFeatureRepresentation *writeFeats = NULL;
  LocalFeatureRepresentation *readFeats = NULL;
  LocalFeatureRepresentation *getFeats = NULL;


  if ( opSiftImpl == "NICE" )
  {
    cSIFT = new LFonHSG ( conf, "HSGtest" );
  }
  else if ( opSiftImpl == "VANDESANDE" )
  {
    // the used features
    cSIFT = new LFColorSande ( conf, "LFColorSandeTest" );
  }
  else
  {
    fthrow ( Exception, "feattype: %s not yet supported" << opSiftImpl );
  }

  getFeats = cSIFT;

  if ( writefeat )
  {
    // write the features to a file, if there isn't any to read
    writeFeats = new LFWriteCache ( conf, cSIFT );
    getFeats = writeFeats;
  }

  if ( readfeat )
  {
    // read the features from a file
    if ( writefeat )
    {
      readFeats = new LFReadCache ( conf, writeFeats, -1 );
    }
    else
    {
      readFeats = new LFReadCache ( conf, cSIFT, -1 );
    }
    getFeats = readFeats;
  }


  // additional Colorfeatures
  LocalFeatureColorWeijer lcw ( conf );

  NICE::ColorImage img;

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

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

  VVector features;
  VVector cfeatures;
  VVector positions;

  getFeats->extractFeatures ( img, features, positions );

  if ( usecolorfeats )
    lcw.getDescriptors ( img, cfeatures, positions );

  set<double> scales;

  int j = 0;
  int lfdimension = -1;
  for ( VVector::const_iterator i = features.begin();
        i != features.end();
        i++, j++ )
  {
    const NICE::Vector & x = *i;

    if ( lfdimension < 0 ) lfdimension = ( int ) x.size();
    else assert ( lfdimension == ( int ) x.size() );

    NICE::Vector *v = new NICE::Vector ( x );

    if ( usecolorfeats )
      v->append ( cfeatures[j] );

    Example tmp = Example ( v );
    tmp.x = ( int ) positions[j][0];
    tmp.y = ( int ) positions[j][1];
    tmp.width = ( int ) ( 16.0 * positions[j][2] );
    tmp.height = tmp.width;
    tmp.scale = positions[j][2];
    scales.insert ( positions[j][2] );
    pce.push_back ( pair<int, Example> ( 0, tmp ) );
  }

  //////////////////
  // PCA anwenden //
  //////////////////
  pce.filename = currentFile;
  if ( usepca )
  {
    doPCA ( pce );
    lfdimension = dim;
  }

  //////////////////
  // BoV anwenden //
  //////////////////
  if ( norm )
    normalize ( pce );
  if ( usegmm || usekmeans )
  {
    if ( !usepca && !norm )
      normalize ( pce );
    convertLowToHigh ( pce );
    smoothHL ( pce );
    lfdimension = gaussians;
  }

  /////////////////////////////////////////
  // Wahrscheinlichkeitskarten erstellen //
  /////////////////////////////////////////
  int klassen = probabilities.channels();
  NICE::MultiChannelImageT<double> preMap ( xsize, ysize, klassen*scales.size() );

  // initialisieren
  for ( int y = 0 ; y < ysize ; y++ )
    for ( int x = 0 ; x < xsize ; x++ )
    {
      // alles zum Hintergrund machen
      segresult.setPixel ( x, y, 0 );
      // Die Wahrscheinlichkeitsmaps auf 0 initialisieren
      for ( int i = 0 ; i < ( int ) probabilities.channels(); i++ )
      {
        probabilities[i](x,y) = 0.0;
      }
      for ( int j = 0; j < ( int ) preMap.channels(); j++ )
      {
        preMap[j](x,y) = 0.0;
      }
    }

  // Die Wahrscheinlichkeitsmaps mit den einzelnen Wahrscheinlichkeiten je Skalierung füllen
  int scalesize = scales.size();

  // Globale Häufigkeiten akkumulieren
  FullVector fV ( ( int ) probabilities.channels() );

  for ( int i = 0; i < fV.size(); i++ )
    fV[i] = 0.0;

  // read allowed classes

  string cndir = conf->gS ( "SemSegCsurka", "cndir", "" );
  int classes = ( int ) probabilities.channels();
  vector<int> useclass ( classes, 1 );

  std::vector< std::string > list;
  StringTools::split ( currentFile, '/', list );

  string orgname = list.back();
  if ( cndir != "" )
  {
    useclass = vector<int> ( classes, 0 );
    ifstream infile ( ( cndir + "/" + orgname + ".dat" ).c_str() );
    while ( !infile.eof() && infile.good() )
    {
      int tmp;
      infile >> tmp;
      if ( tmp >= 0 && tmp < classes )
      {
        useclass[tmp] = 1;
      }
    }
  }

#ifdef UNCERTAINTY
  std::vector<FloatImage> uncert;
  std::vector<FloatImage> gpUncertainty;
  std::vector<FloatImage> gpMean;    
  std::vector<FloatImage> gpMeanRatio;  
  std::vector<FloatImage> gpWeightAll;
  std::vector<FloatImage> gpWeightRatio;
//   std::vector<FloatImage> gpImpactAll;
//   std::vector<FloatImage> gpImpactRatio;
  
  //pre-allocate storage -- one image per scale and method
  for(int s = 0; s < scalesize; s++)
  {
    uncert.push_back(FloatImage(xsize, ysize));
    uncert[s].set(0.0);
    
    gpUncertainty.push_back(FloatImage(xsize, ysize));
    gpMean.push_back(FloatImage(xsize, ysize));
    gpMeanRatio.push_back(FloatImage(xsize, ysize));
    gpWeightAll.push_back(FloatImage(xsize, ysize));
    gpWeightRatio.push_back(FloatImage(xsize, ysize));
/*    gpImpactAll.push_back(FloatImage(xsize, ysize));    
    gpImpactRatio.push_back(FloatImage(xsize, ysize));   */  
   
    gpUncertainty[s].set(0.0);
    gpMean[s].set(0.0);
    gpMeanRatio[s].set(0.0);
    gpWeightAll[s].set(0.0);
    gpWeightRatio[s].set(0.0);
//     gpImpactAll[s].set(0.0); 
//     gpImpactRatio[s].set(0.0);   
  }
  
  ColorImage imgrgb ( xsize, ysize );
  std::string s;
  std::stringstream out;
  std::vector< std::string > list2;
  StringTools::split ( Globals::getCurrentImgFN (), '/', list2 );
  out << "uncertainty/" << list2.back();
  
  double maxu = -numeric_limits<double>::max();
  double minu = numeric_limits<double>::max();
  
  double gpNoise =  conf->gD("GPHIK", "noise", 0.01);
  
#endif

  #ifdef UNCERTAINTY
  std::cerr << "compute values for uncertainty stuff as well" << std::endl;
  #endif
  
  if ( classifier != NULL )
  {
    clog << "[log] SemSegCsruka::classifyregions: Wahrscheinlichkeitskarten erstellen: classifier != NULL" << endl;
#pragma omp parallel for
    for ( int s = 0; s < scalesize; s++ )
    {
#pragma omp parallel for
      for ( int i = s; i < ( int ) pce.size(); i += scalesize )
      {
        ClassificationResult r = classifier->classify ( pce[i].second );

        #ifdef UNCERTAINTY
        //we need this if we want to compute GP-AL-measure lateron
        double minMeanAbs ( numeric_limits<double>::max() );
        double maxMeanAbs ( 0.0 );
        double sndMaxMeanAbs ( 0.0 );
        #endif
        
        for ( int j = 0 ; j < r.scores.size(); j++ )
        {
          if ( useclass[j] == 0 )
            continue;

          fV[j] += r.scores[j];
          preMap.set ( pce[i].second.x, pce[i].second.y, r.scores[j], j + s*klassen );
          
         #ifdef UNCERTAINTY 
          //check whether we found a class with higher smaller abs mean than the current minimum
         if (abs(r.scores[j]) < minMeanAbs)  
           minMeanAbs = abs(r.scores[j]);
         //check for larger abs mean as well
         if (abs(r.scores[j]) > maxMeanAbs)
         {
           sndMaxMeanAbs = maxMeanAbs;
           maxMeanAbs = abs(r.scores[j]);
         }
         // and also for the second highest mean of all classes
         else if (abs(r.scores[j]) > sndMaxMeanAbs)
         {
           sndMaxMeanAbs = abs(r.scores[j]);
         }
         #endif          
        }

        /*if(r.uncertainty < 0.0)
        {
          cerr << "uncertainty: " << r.uncertainty << endl;
          pce[i].second.svec->store(cerr);
          cerr << endl;
          exit(-1);
        }*/
#ifdef UNCERTAINTY
        uncert[s] ( pce[i].second.x, pce[i].second.y ) = r.uncertainty;
        maxu = std::max ( r.uncertainty, maxu );
        minu = std::min ( r.uncertainty, minu );
        
        
        double firstTerm (1.0 / sqrt(r.uncertainty+gpNoise));
        
        //compute the heuristic GP-UNCERTAINTY, as proposed by Kapoor et al. in IJCV 2010
        // GP-UNCERTAINTY : |mean| / sqrt(var^2 + gpnoise^2)
        gpUncertainty[s] ( pce[i].second.x, pce[i].second.y ) = maxMeanAbs*firstTerm; //firstTerm = 1.0 / sqrt(r.uncertainty+gpNoise))
        
        // compute results when we take the lowest mean value of all classes
        gpMean[s] ( pce[i].second.x, pce[i].second.y ) = minMeanAbs;
        
        //look at the difference in the absolut mean values for the most plausible class
        // and the second most plausible class
        gpMeanRatio[s] ( pce[i].second.x, pce[i].second.y ) = maxMeanAbs - sndMaxMeanAbs;
        

        //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 < r.scores.size(); j++ )
        {
          if ( useclass[j] == 0 )
            continue;
          // look at the difference to plus 1          
          diffToPositive.push_back(abs(r.scores[j] - 1));
          // look at the difference to -1          
          diffToNegative.push_back(abs(r.scores[j] + 1));
          //sum up the difference to -1
          diffToNegativeSum += abs(r.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;
        }
        gpWeightAll[s] ( pce[i].second.x, pce[i].second.y ) = thirdTerm*firstTerm;        
        
        //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
        gpWeightRatio[s] ( pce[i].second.x, pce[i].second.y ) = (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;      
#endif
      }
    }
  }
  else
  {
//#pragma omp parallel for
    for ( int s = 0; s < scalesize; s++ )
    {
//#pragma omp parallel for
      for ( int i = s; i < ( int ) pce.size(); i += scalesize )
      {
        ClassificationResult r = vclassifier->classify ( * ( pce[i].second.vec ) );
        
        #ifdef UNCERTAINTY
        //we need this if we want to compute GP-AL-measure lateron
        double minMeanAbs ( numeric_limits<double>::max() );
        double maxMeanAbs ( 0.0 );
        double sndMaxMeanAbs ( 0.0 );
        #endif        
        
        for ( int j = 0 ; j < ( int ) r.scores.size(); j++ )
        {
          if ( useclass[j] == 0 )
            continue;
          fV[j] += r.scores[j];
          preMap.set ( pce[i].second.x, pce[i].second.y, r.scores[j], j + s*klassen );
          
         #ifdef UNCERTAINTY 
          //check whether we found a class with higher smaller abs mean than the current minimum
         if (abs(r.scores[j]) < minMeanAbs)  
           minMeanAbs = abs(r.scores[j]);
         //check for larger abs mean as well
         if (abs(r.scores[j]) > maxMeanAbs)
         {
           sndMaxMeanAbs = maxMeanAbs;
           maxMeanAbs = abs(r.scores[j]);
         }
         // and also for the second highest mean of all classes
         else if (abs(r.scores[j]) > sndMaxMeanAbs)
         {
           sndMaxMeanAbs = abs(r.scores[j]);
         }
         #endif            
        }
#ifdef UNCERTAINTY
        uncert[s] ( pce[i].second.x, pce[i].second.y ) = r.uncertainty;
        maxu = std::max ( r.uncertainty, maxu );
        minu = std::min ( r.uncertainty, minu );
        
        
        double firstTerm (1.0 / sqrt(r.uncertainty+gpNoise));
        
        //compute the heuristic GP-UNCERTAINTY, as proposed by Kapoor et al. in IJCV 2010
        // GP-UNCERTAINTY : |mean| / sqrt(var^2 + gpnoise^2)
        gpUncertainty[s] ( pce[i].second.x, pce[i].second.y ) = maxMeanAbs*firstTerm; //firstTerm = 1.0 / sqrt(r.uncertainty+gpNoise))
        
        // compute results when we take the lowest mean value of all classes
        gpMean[s] ( pce[i].second.x, pce[i].second.y ) = minMeanAbs;
        
        //look at the difference in the absolut mean values for the most plausible class
        // and the second most plausible class
        gpMeanRatio[s] ( pce[i].second.x, pce[i].second.y ) = maxMeanAbs - sndMaxMeanAbs;
        

        //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 < fV.size(); j++ )
        {
          if ( useclass[j] == 0 )
            continue;
          // look at the difference to plus 1          
          diffToPositive.push_back(abs(r.scores[j] - 1));
          // look at the difference to -1          
          diffToNegative.push_back(abs(r.scores[j] + 1));
          //sum up the difference to -1
          diffToNegativeSum += abs(r.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;
        }
        gpWeightAll[s] ( pce[i].second.x, pce[i].second.y ) = thirdTerm*firstTerm;        
        
        //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
        gpWeightRatio[s] ( pce[i].second.x, pce[i].second.y ) = (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;      
#endif
      }
    }
  }

  #ifdef UNCERTAINTY
  std::cerr << "uncertainty values and derived scores successfully computed" << std::endl;
  #endif

#ifdef UNCERTAINTY
  cout << "maxvdirect: " << maxu << " minvdirect: " << minu << endl;
  //pre-allocate the image for filtering lateron
  FloatImage gaussUncert ( xsize, ysize );
  
  //just store the first scale
  ICETools::convertToRGB ( uncert[0], imgrgb );
  imgrgb.write ( out.str() + "rough.ppm" );
  
  //pre-allocate memory for filtering of scales
  FloatImage gaussGPUncertainty ( xsize, ysize );
  FloatImage gaussGPMean ( xsize, ysize );
  FloatImage gaussGPMeanRatio( xsize, ysize );
  FloatImage gaussGPWeightAll ( xsize, ysize );
  FloatImage gaussGPWeightRatio ( xsize, ysize );
   
  //just store the first scale for every method
  ICETools::convertToRGB ( gpUncertainty[0], imgrgb );
  imgrgb.write ( out.str() + "gpUncertainty.ppm" );
  ICETools::convertToRGB ( gpMean[0], imgrgb );
  imgrgb.write ( out.str() + "gpMean.ppm" );
  ICETools::convertToRGB ( gpMeanRatio[0], imgrgb );
  imgrgb.write ( out.str() + "gpMeanRatio.ppm" );
  ICETools::convertToRGB ( gpWeightAll[0], imgrgb );
  imgrgb.write ( out.str() + "gpWeightAll.ppm" );
  ICETools::convertToRGB ( gpWeightRatio[0], imgrgb );
  imgrgb.write ( out.str() + "gpWeightRatio.ppm" );  
  
#endif

  vector<double> scalesVec;
  for ( set<double>::const_iterator iter = scales.begin();
        iter != scales.end();
        ++iter )
  {
    scalesVec.push_back ( *iter );
  }

#undef VISSEMSEG
#ifdef VISSEMSEG

  for ( int j = 0 ; j < ( int ) preMap.channels(); j++ )
  {
    cout << "klasse: " << j << endl;//" " << cn.text ( j ) << endl;

    NICE::Matrix tmp ( preMap.ysize, preMap.xsize );
    double maxval = 0.0;
    for ( int y = 0; y < preMap.ysize; y++ )
      for ( int x = 0; x < preMap.xsize; x++ )
      {
        double val = preMap.get ( x, y, j );
        tmp ( y, x ) = val;
        maxval = std::max ( val, maxval );
      }

    NICE::ColorImage imgrgb ( preMap.xsize, preMap.ysize );
    ICETools::convertToRGB ( tmp, imgrgb );

    cout << "maxval = " << maxval << " for class " << j << endl; //cn.text ( j ) << endl;

    //Show ( ON, imgrgb, cn.text ( j ) );
    //showImage(imgrgb, "Ergebnis");

    std::string s;
    std::stringstream out;
    out << "tmpprebmap" << j << ".ppm";
    s = out.str();
    imgrgb.writePPM ( s );

    //getchar();
  }
#endif

  // Gaußfiltern
  clog << "[log] SemSegCsruka::classifyregions: Wahrscheinlichkeitskarten erstellen -> Gaussfiltern" << endl;
  for ( int s = 0; s < scalesize; s++ )
  {
    double sigma = sigmaweight * 16.0 * scalesVec[s];
    cerr << "sigma: " << sigma << endl;
#pragma omp parallel for
    for ( int i = 0; i < klassen; i++ )
    {
      if ( forbidden_classes.find ( i ) != forbidden_classes.end() )
      {
        continue;
      }

      int pos = i + s * klassen;

      double maxval = preMap[pos](0,0);
      double minval = maxval;

      for ( int y = 0; y < ysize; y++ )
      {
        for ( int x = 0; x < xsize; x++ )
        {
          maxval = std::max ( maxval, preMap[pos](x,y) );
          minval = std::min ( minval, preMap[pos](x,y) );
        }
      }

      NICE::FloatImage dblImg ( xsize, ysize );
      NICE::FloatImage gaussImg ( xsize, ysize );

      for ( int y = 0; y < ysize; y++ )
      {
        for ( int x = 0; x < xsize; x++ )
        {
          dblImg.setPixel ( x, y, preMap[pos](x,y) );
        }
      }

      NICE::FilterT<float, float, float> filt;
      filt.filterGaussSigmaApproximate ( dblImg, sigma, &gaussImg );

      for ( int y = 0; y < ysize; y++ )
      {
        for ( int x = 0; x < xsize; x++ )
        {
          preMap[pos](x,y) = gaussImg.getPixel ( x, y );
        }
      }
    }
#ifdef UNCERTAINTY
    filterGaussSigmaApproximate<float, float, float> ( uncert[s], sigma, &gaussUncert );
    uncert[s] = gaussUncert;
    
    //apply the gauss-filtering to all scales of every method
    filterGaussSigmaApproximate<float, float, float> ( gpUncertainty[s], sigma, &gaussGPUncertainty );
    filterGaussSigmaApproximate<float, float, float> ( gpMean[s], sigma, &gaussGPMean );
    filterGaussSigmaApproximate<float, float, float> ( gpMeanRatio[s], sigma, &gaussGPMeanRatio );
    filterGaussSigmaApproximate<float, float, float> ( gpWeightAll[s], sigma, &gaussGPWeightAll );
    filterGaussSigmaApproximate<float, float, float> ( gpWeightRatio[s], sigma, &gaussGPWeightRatio );
    
    gpUncertainty[s] = gaussGPUncertainty; 
    gpMean[s] = gaussGPMean; 
    gpMeanRatio[s] = gaussGPMeanRatio; 
    gpWeightAll[s] = gaussGPWeightAll;
    gpWeightRatio[s] = gaussGPWeightRatio;   
#endif
  }

  // Zusammenfassen und auswerten
  clog << "[log] SemSegCsruka::classifyregions: Wahrscheinlichkeitskarten erstellen -> zusammenfassen und auswerten" << endl;
//#pragma omp parallel for
  for ( int x = 0; x < xsize; x++ )
  {
    for ( int y = 0; y < ysize; y++ )
    {
      for ( int j = 0 ; j < ( int ) probabilities.channels(); j++ )
      {
        double prob = 0.0;
        for ( int s = 0; s < ( int ) scalesize; s++ )
        {
          prob += preMap.get ( x, y, j + s * klassen );
        }

        double val = prob / ( double ) ( scalesize );
        probabilities.set ( x, y, val, j );
      }
    }
  }
  
#ifdef UNCERTAINTY
  for ( int x = 0; x < xsize; x++ )
  {
    for ( int y = 0; y < ysize; y++ )
    {
      for ( int s = 0; s < ( int ) scalesize; s++ )
      {
        gaussUncert(x,y) += uncert[s](x,y);
        //and for the other methods as well
        gaussGPUncertainty(x,y) += gpUncertainty[s](x,y);
        gaussGPMean(x,y) += gpMean[s](x,y);
        gaussGPMeanRatio(x,y) += gpMeanRatio[s](x,y);
        gaussGPWeightAll(x,y) += gpWeightAll[s](x,y);
        gaussGPWeightRatio(x,y) += gpWeightRatio[s](x,y);
      }
      gaussUncert(x,y)/=scalesize;
      //and for the other methods as well
      gaussGPUncertainty(x,y)/=scalesize;
      gaussGPMean(x,y)/=scalesize;
      gaussGPMeanRatio(x,y)/=scalesize;
      gaussGPWeightAll(x,y)/=scalesize;
      gaussGPWeightRatio(x,y)/=scalesize;      
    }
  }

  maxu = -numeric_limits<double>::max();
  minu = numeric_limits<double>::max();
  for ( int y = 0; y < ysize; y++ )
  {
    for ( int x = 0; x < xsize; x++ )
    {
      double val = uncert[0] ( x, y );
      maxu = std::max ( val, maxu );
      minu = std::min ( val, minu );
    }
  }
  cout << "maxvo = " << maxu << " minvo = " << minu << endl;

  maxu = -numeric_limits<float>::max();
  minu = numeric_limits<float>::max();

  for ( int y = 0; y < ysize; y++ )
  {
    for ( int x = 0; x < xsize; x++ )
    {
      double val = gaussUncert ( x, y );
      maxu = std::max ( val, maxu );
      minu = std::min ( val, minu );
    }
  }
  cout << "maxvf = " << maxu << " minvf = " << minu << endl;
  
  gaussUncert(0,0) = 0.0;
  gaussUncert(0,1) = 0.04;
  ICETools::convertToRGB ( gaussUncert, imgrgb );
  imgrgb.write ( out.str() + "filtered.ppm" );
  
  ICETools::convertToRGB ( gaussGPUncertainty, imgrgb );
  imgrgb.write ( out.str() + "gpUncertaintyFiltered.ppm" );
  ICETools::convertToRGB ( gaussGPMean, imgrgb );
  imgrgb.write ( out.str() + "gpMeanFiltered.ppm" );
  ICETools::convertToRGB ( gaussGPMeanRatio, imgrgb );
  imgrgb.write ( out.str() + "gpMeanRatioFiltered.ppm" );
  ICETools::convertToRGB ( gaussGPWeightAll, imgrgb );
  imgrgb.write ( out.str() + "gpWeightAllFiltered.ppm" );
  ICETools::convertToRGB ( gaussGPWeightRatio, imgrgb );
  imgrgb.write ( out.str() + "gpWeightRatioFiltered.ppm" );  
  
#endif

#undef VISSEMSEG
#ifdef VISSEMSEG

  std::string s;
  std::stringstream out;
  std::vector< std::string > list2;
  StringTools::split ( Globals::getCurrentImgFN (), '/', list2 );

  out << "probmaps/" << list2.back() << ".probs";

  s = out.str();

  probabilities.store ( s );

  for ( int j = 0 ; j < ( int ) probabilities.channels(); j++ )
  {
    cout << "klasse: " << j << endl;//" " << cn.text ( j ) << endl;

    NICE::Matrix tmp ( probabilities.ysize, probabilities.xsize );
    double maxval = 0.0;
    for ( int y = 0; y < probabilities.ysize; y++ )
      for ( int x = 0; x < probabilities.xsize; x++ )
      {
        double val = probabilities.get ( x, y, j );

        tmp ( y, x ) = val;
        maxval = std::max ( val, maxval );
      }

    NICE::ColorImage imgrgb ( probabilities.xsize, probabilities.ysize );
    ICETools::convertToRGB ( tmp, imgrgb );

    cout << "maxval = " << maxval << " for class " << j << endl; //cn.text ( j ) << endl;

    //Show ( ON, imgrgb, cn.text ( j ) );
    //showImage(imgrgb, "Ergebnis");

    std::string s;
    std::stringstream out;
    out << "tmp" << j << ".ppm";
    s = out.str();
    imgrgb.writePPM ( s );

    //getchar();
  }
#endif
  if ( useregions )
  {
    if ( bestclasses > 0 )
    {
      PSSImageLevelPrior pss ( 0, bestclasses, 0.2 );
      pss.setPrior ( fV );
      pss.postprocess ( segresult, probabilities );
    }

    //Regionen ermitteln

    int regionsize = seg->segRegions ( img, mask );

    Regionen.clear();
    vector<vector <double> > regionprob;

#ifdef UNCERTAINTY
    std::vector<double> regionUncert;
    
    std::vector<double> regionGPUncertainty;
    std::vector<double> regionGPMean;
    std::vector<double> regionGPMeanRatio;
    std::vector<double> regionGPWeightAll;
    std::vector<double> regionGPWeightRatio;    
#endif

    // Wahrscheinlichkeiten für jede Region initialisieren
    for ( int i = 0; i < regionsize; i++ )
    {
      vector<double> tmp;
      for ( int j = 0; j < ( int ) probabilities.channels(); j++ )
      {
        tmp.push_back ( 0.0 );
      }
      regionprob.push_back ( tmp );
      Regionen.push_back ( pair<int, Example> ( 0, Example() ) );
#ifdef UNCERTAINTY
      regionUncert.push_back ( 0.0 );
      
      regionGPUncertainty.push_back ( 0.0 );
      regionGPMean.push_back ( 0.0 );
      regionGPMeanRatio.push_back ( 0.0 );
      regionGPWeightAll.push_back ( 0.0 );
      regionGPWeightRatio.push_back ( 0.0 );
#endif
    }

    // Wahrscheinlichkeiten für Regionen bestimmen
    for ( int x = 0; x < xsize; x++ )
    {
      for ( int y = 0; y < ysize; y++ )
      {
        int pos = mask ( x, y );
        Regionen[pos].second.weight += 1.0;
        Regionen[pos].second.x += x;
        Regionen[pos].second.y += y;
        for ( int j = 0 ; j < ( int ) probabilities.channels(); j++ )
        {
          double val = probabilities.get ( x, y, j );
          regionprob[pos][j] += val;
        }
#ifdef UNCERTAINTY
        regionUncert[pos] += gaussUncert ( x, y );
        
        regionGPUncertainty[pos] += gaussGPUncertainty ( x, y );
        regionGPMean[pos] += gaussGPMean ( x, y );
        regionGPMeanRatio[pos] += gaussGPMeanRatio ( x, y );
        regionGPWeightAll[pos] += gaussGPWeightAll ( x, y );
        regionGPWeightRatio[pos] += gaussGPWeightRatio ( x, y );
#endif
      }
    }

    /*
    cout << "regions: " << regionsize << endl;
    cout << "outfeats: " << endl;
    for(int j = 0; j < regionprob.size(); j++)
    {
     for(int i = 0; i < regionprob[j].size(); i++)
     {
      cout << regionprob[j][i] << " ";
     }
     cout << endl;
    }
    cout << endl;
    getchar();*/

    // beste Wahrscheinlichkeit je Region wählen
    for ( int i = 0; i < regionsize; i++ )
    {
      if ( Regionen[i].second.weight > 0 )
      {
        Regionen[i].second.x /= ( int ) Regionen[i].second.weight;
        Regionen[i].second.y /= ( int ) Regionen[i].second.weight;
      }
      double maxval = -numeric_limits<double>::max();
      int maxpos = 0;

      for ( int j = 0 ; j < ( int ) regionprob[i].size(); j++ )
      {
        if ( forbidden_classes.find ( j ) != forbidden_classes.end() )
          continue;

        regionprob[i][j] /= Regionen[i].second.weight;

        if ( maxval < regionprob[i][j] )
        {
          maxval = regionprob[i][j];
          maxpos = j;
        }
        probabilities.set ( Regionen[i].second.x, Regionen[i].second.y, regionprob[i][j], j );
      }

      Regionen[i].first = maxpos;
#ifdef UNCERTAINTY
      regionUncert[i] /= Regionen[i].second.weight;
      
      regionGPUncertainty[i] /= Regionen[i].second.weight;
      regionGPMean[i] /= Regionen[i].second.weight;
      regionGPMeanRatio[i] /= Regionen[i].second.weight;
      regionGPWeightAll[i] /= Regionen[i].second.weight;
      regionGPWeightRatio[i] /= Regionen[i].second.weight;
#endif
    }
    // Pixel jeder Region labeln
    for ( int y = 0; y < ( int ) mask.cols(); y++ )
    {
      for ( int x = 0; x < ( int ) mask.rows(); x++ )
      {
        int pos = mask ( x, y );
        segresult.setPixel ( x, y, Regionen[pos].first );
#ifdef UNCERTAINTY
        gaussUncert ( x, y ) = regionUncert[pos];
        
        gaussGPUncertainty ( x, y ) = regionGPUncertainty[pos];
        gaussGPMean ( x, y ) = regionGPMean[pos];
        gaussGPMeanRatio ( x, y ) = regionGPMeanRatio[pos];
        gaussGPWeightAll ( x, y ) = regionGPWeightAll[pos];
        gaussGPWeightRatio ( x, y ) = regionGPWeightRatio[pos];        
#endif
      }
    }   
#ifdef UNCERTAINTY
    maxu = -numeric_limits<float>::max();
    minu = numeric_limits<float>::max();
    for ( int y = 0; y < ysize; y++ )
    {
      for ( int x = 0; x < xsize; x++ )
      {
        //float val = uncert(x,y);
        double val = gaussUncert ( x, y );
        maxu = std::max ( val, maxu );
        minu = std::min ( val, minu );
      }
    }
    cout << "maxvr = " << maxu << " minvr = " << minu << endl;
//    uncert(0,0) = 1;
//    uncert(0,1) = 0;
    ICETools::convertToRGB ( gaussUncert, imgrgb );
    imgrgb.write ( out.str() + "region.ppm" );
    
  ICETools::convertToRGB ( gaussGPUncertainty, imgrgb );
  imgrgb.write ( out.str() + "gpUncertaintyRegion.ppm" );
  ICETools::convertToRGB ( gaussGPMean, imgrgb );
  imgrgb.write ( out.str() + "gpMeanRegion.ppm" );
  ICETools::convertToRGB ( gaussGPMeanRatio, imgrgb );
  imgrgb.write ( out.str() + "gpMeanRatioRegion.ppm" );
  ICETools::convertToRGB ( gaussGPWeightAll, imgrgb );
  imgrgb.write ( out.str() + "gpWeightAllRegion.ppm" );
  ICETools::convertToRGB ( gaussGPWeightRatio, imgrgb );
  imgrgb.write ( out.str() + "gpWeightRatioRegion.ppm" );      
#endif

#undef WRITEREGIONS
#ifdef WRITEREGIONS
    RegionGraph rg;
    seg->getGraphRepresentation ( img, mask, rg );
    for ( uint pos = 0; pos < regionprob.size(); pos++ )
    {
      rg[pos]->setProbs ( regionprob[pos] );
    }

    std::string s;
    std::stringstream out;
    std::vector< std::string > list;
    StringTools::split ( Globals::getCurrentImgFN (), '/', list );

    out << "rgout/" << list.back() << ".graph";
    string writefile = out.str();
    rg.write ( writefile );
#endif
  }
  else
  {

    PSSImageLevelPrior pss ( 1, 4, 0.2 );
    pss.setPrior ( fV );
    pss.postprocess ( segresult, probabilities );

  }

  // Saubermachen:
  clog << "[log] SemSegCsurka::classifyregions: sauber machen" << endl;
  for ( int i = 0; i < ( int ) pce.size(); i++ )
  {
    pce[i].second.clean();
  }
  pce.clear();

  if ( cSIFT != NULL )
    delete cSIFT;
  if ( writeFeats != NULL )
    delete writeFeats;
  if ( readFeats != NULL )
    delete readFeats;
  getFeats = NULL;
}

void SemSegCsurka::semanticseg (
        CachedExample *ce,
        NICE::ImageT<int> & segresult,
        NICE::MultiChannelImageT<double> & probabilities )
{
  Examples regions;
  NICE::Matrix regionmask;
  classifyregions ( ce, segresult, probabilities, regions, regionmask );
  if ( userellocprior || srg != NULL || gcopt != NULL )
  {
    if ( userellocprior )
      relloc->postprocess ( regions, probabilities );

    if ( srg != NULL )
      srg->optimizeShape ( regions, regionmask, probabilities );

    if ( gcopt != NULL )
      gcopt->optimizeImage ( regions, regionmask, probabilities );

    // Pixel jeder Region labeln
    for ( int y = 0; y < ( int ) regionmask.cols(); y++ )
    {
      for ( int x = 0; x < ( int ) regionmask.rows(); x++ )
      {
        int pos = regionmask ( x, y );
        segresult.setPixel ( x, y, regions[pos].first );
      }
    }
  }

#ifndef NOVISUAL
#undef VISSEMSEG
#ifdef VISSEMSEG
//  showImage(img);
  for ( int j = 0 ; j < ( int ) probabilities.channels(); j++ )
  {
    cout << "klasse: " << j << " " << cn.text ( j ) << endl;

    NICE::Matrix tmp ( probabilities.ysize, probabilities.xsize );
    double maxval = -numeric_limits<double>::max();
    for ( int y = 0; y < probabilities.ysize; y++ )
      for ( int x = 0; x < probabilities.xsize; x++ )
      {
        double val = probabilities.get ( x, y, j );
        tmp ( y, x ) = val;
        maxval = std::max ( val, maxval );
      }

    NICE::ColorImage imgrgb ( probabilities.xsize, probabilities.ysize );
    ICETools::convertToRGB ( tmp, imgrgb );

    cout << "maxval = " << maxval << " for class " << cn.text ( j ) << endl;

    Show ( ON, imgrgb, cn.text ( j ) );
    imgrgb.Write ( "tmp.ppm" );

    getchar();
  }
#endif
#endif

}