#include <SizeCluster_SVM.hpp>
//#include <Progress.hpp>
#include <CellCluster.hpp>
#include <math.h> 
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_vector.h>
#include <gsl/gsl_eigen.h>
#include <dlib/statistics.h>
#include <iostream>
#include <vector>
#include <fstream>
#include <dlib/mlp.h>
#include <limits>
#include <Progress.hpp>
#include <QFileDialog>
#include <QFileDialog>
#include <dlib/svm.h>
using namespace dlib;

namespace mgx
{
  
  std::map<Point3d,Point3d> MapPosAvg; 
  int label, nhbrlabel;
  typedef std::pair<int, double> loopMaps;
  typedef std::pair<VertexPr, Set> VertexPrSet;
  typedef std::pair<VertexPr, Point3d> VertexPrPoint;
  std::map<int,int> countTempo, Totneighbors;
  std::map<vertex,int>countavg, VContain;
  #define PI M_PI
  /*Paramter calculations for each selected cell*/
  //Calculate Area
  void CalculateArea(vvGraph &T, CellCluster &cc){  
    forall(const vertex &v, T){  
      if(v->label == -1)
        continue;
      forall(const vertex &n, T.neighbors(v)){
        vertex m = T.nextTo(v,n);
        cc.Area[v->label] += triangleArea(v->pos, n->pos, m->pos);     
      }
      std::cout<<"\nArea\t"<<cc.Area[v->label];
    }
    cc.FeatureVector.push_back(cc.Area);
  }
  //Calculate Perimeter
  void CalculatePerimeter(vvGraph &T, CellCluster &cc){
    forall(const vertex &v, T){
      if(v->label == -1)
        continue;
      forall(const vertex &n, T.neighbors(v)){
        vertex m = T.nextTo(v, n);
        if(T.edge(m,n) or T.edge(m, n))
          cc.Perimeter[v->label] += norm(n->pos-m->pos);
      }
    }
    cc.FeatureVector.push_back(cc.Perimeter);
  }
  //Calculate bending - the sum of cross products of adjecent edges in a cell 
  void CalculateBending(vvGraph &T, CellCluster &cc){ 
    forall(const vertex &v, T){
      if(v->label == -1)
        continue;
      forall(const vertex &n, T.neighbors(v)){
        vertex m = T.nextTo(v, n);
        vertex l = T.nextTo(v, m);
        Point3d nm =  m->pos - n->pos;
        nm /= norm(nm);
        Point3d ml =  l->pos - m->pos;
        ml /= norm(ml);
        float s = norm(nm ^ ml);
        float c = (nm * ml);
        if(c < 0)
          s = 1 + (1-s);
        cc.Bending[v->label] += s;
      }
    }
    cc.FeatureVector.push_back(cc.Bending);
  }
  //Calculate neighborhood bending using closest neighbor map
  void CalculateStomatalBending(vvGraph &T, CellCluster &cc){ 
    forall(VertexPr Pr, cc.NeighbourMap) {  
      double stomataBending = 0;
      if(Pr.first->label == 27 or Pr.first->label == 55)
        std::cout<<"\n"<<Pr.first->label<<"\t"<<Pr.second->label <<"t"<<cc.NhbrElem[Pr].size();
      if(cc.NhbrElem[Pr].size() <= 1)
        continue;
      forall(const vertex &n, T.neighbors(Pr.first)){
        vertex m = T.nextTo(Pr.first, n);
        vertex l = T.nextTo(Pr.first, m);
        if(cc.NhbrElem[Pr].find(m) != cc.NhbrElem[Pr].end() and cc.NhbrElem[Pr].find(l) != cc.NhbrElem[Pr].end()){   
          Point3d nm =  m->pos - n->pos;
          nm /= norm(nm);
          Point3d ml =  l->pos - m->pos;
          ml /= norm(ml);
          float s = norm(nm ^ ml);
          float c = (nm * ml);
          if(c < 0)
            s = 1 + (1-s);
          stomataBending += s;
        }
      }  
      forall(const vertex &n, T.neighbors(Pr.second)){
        vertex m = T.nextTo(Pr.second, n);
        vertex l = T.nextTo(Pr.second, m);
        if(cc.NhbrElem[Pr].find(m) != cc.NhbrElem[Pr].end() and cc.NhbrElem[Pr].find(l) != cc.NhbrElem[Pr].end()){   
          Point3d nm =  m->pos - n->pos;
          nm /= norm(nm);
          Point3d ml =  l->pos - m->pos;
          ml /= norm(ml);
          float s = norm(nm ^ ml);
          float c = (nm * ml);
          if(c < 0)
            s = 1 + (1-s);
          stomataBending += s;
        }
      }
      cc.stomataBending[Pr.first->label] = cc.Bending[Pr.first->label] + cc.Bending[Pr.second->label] - stomataBending;
      cc.stomataBending[Pr.second->label] = cc.Bending[Pr.first->label] + cc.Bending[Pr.second->label] - stomataBending; 
      //std::cout<<"\nfor\t"<<Pr.first->label<<"\t"<<Pr.second->label<<"\tstomatal bending\t"<<cc.stomataBending[Pr.first->label]<<"\t"<<cc.stomataBending[Pr.second->label]<<"\t"<<stomataBending<<"\t"<<cc.NhbrElem[Pr].size();
    }
    cc.FeatureVector.push_back(cc.stomataBending);
  }
  void CalculateLenBr(vvGraph &T, CellCluster &cc){
    forall(const vertex &v, T){
      if(v->label == -1)
        continue;
      float min = 0;
      float max = 1000;
      forall(const vertex &n, T.neighbors(v)){
        vertex m = T.nextTo(v, n); 
        if((norm(v->pos - n->pos)) > min){
          min = norm(v->pos - n->pos); 
          cc.ClusterLength[v->label] = min;  
        }
        if((norm(v->pos - n->pos)) < max){
          max = norm(v->pos - n->pos); 
          cc.ClusterBreadth[v->label] = max;  
        } 
      }
      cc.DiffLenBr[v->label] = cc.ClusterLength[v->label]/cc.ClusterBreadth[v->label];
    }
    cc.FeatureVector.push_back(cc.DiffLenBr);            
  }
  void CalculateNeighbours(vvGraph &T, CellCluster &cc){  
    forall(const vertex &v, T){
      if(v->label == -1)
        continue;
      forall(const vertex &n, T.neighbors(v)){
        if(n->label == -1){
          forall(const vertex &q, T.neighbors(n)){
            if(q->label != -1 and q->label != v->label and cc.Nhbr[v->label].find(q) == cc.Nhbr[v->label].end())
              cc.Nhbr[v->label].insert(q);
          }
        }
      }
      cc.TotalNeighbours[v->label] = cc.Nhbr[v->label].size();
    } 
    cc.FeatureVector.push_back(cc.TotalNeighbours);
  }
  void CalculateCurvature(vvGraph &T, CellCluster &cc){  
    forall(const vertex &v, T){
      if(v->label == -1)
        continue;
      cc.Curvature[v->label] = pow(cc.Perimeter[v->label], 2)/cc.Area[v->label];
    }
    cc.FeatureVector.push_back(cc.Curvature);
  }
  void FindCommonNeighbours(vvGraph &T, CellCluster &cc){
    std::map<VertexPr, Set> NhbrElem;
    forall(const vertex &v, T) {
      if(v->label == -1)
        continue;
      forall(const vertex &n, cc.Nhbr[v->label]){
        if(n->label == -1)
          continue;
        if(NhbrElem.find(std::make_pair(v, n)) != NhbrElem.end() or NhbrElem.find(std::make_pair(n, v)) != NhbrElem.end())
          continue;
        if(v->label == 2 and n->label == 3)
          std::cout<<"\ncond\t"<<(NhbrElem.find(std::make_pair(v, n)) != NhbrElem.end())<<"\t"<<(NhbrElem.find(std::make_pair(n, v)) != NhbrElem.end());       
        forall(const vertex &x, T.neighbors(v)){  //Iterate through the neighbours of v
          forall(const vertex &y, T.neighbors(n)){  //Iterate through the neighbours of n
            if(x == y and x != n and y != v){           
              NhbrElem[std::make_pair(v, n)].insert(x);
              NhbrElem[std::make_pair(v, n)].insert(y);
            }
          }
        }
        if(NhbrElem[std::make_pair(v, n)].size() > 1){
          cc.NhbrElem[std::make_pair(v, n)] = NhbrElem[std::make_pair(v, n)];
          std::cout<<"\nSize\t"<<v->label<<"\t"<<n->label<<"\t"<<cc.NhbrElem[std::make_pair(v, n)].size()<<"\t"<<NhbrElem[std::make_pair(v, n)].size();
        }
      } 
    }
    std::cout<<"\ncc.NhbrElem\t"<<cc.NhbrElem.size();  
  }
  void CalculateStandardDeviation(vvGraph &T, CellCluster &cc){  
    forall(VertexPrSet CNhbr, cc.NhbrElem){  //Calculate center of the border of neighbouring vertices
      cc.meanNhbr[CNhbr.first] = 0;  //cc.NhbrElem contains common neighbours (border vertices)
      forall(const vertex &elem , CNhbr.second){  //Iterate through neighbours
        cc.meanNhbr[CNhbr.first] += elem->pos;  //Add positions
      } 
      cc.meanNhbr[CNhbr.first] /= CNhbr.second.size();  //Take average of positions
      cc.StDev[CNhbr.first.first->label] = HUGE_VAL;  
      cc.StDev[CNhbr.first.second->label] = HUGE_VAL;
    }
    forall(VertexPrSet CNhbr, cc.NhbrElem){  //
      double stdev = 0, countele = 0;  
        forall(const vertex &elem, T.neighbors(CNhbr.first.first)){
          if(CNhbr.second.find(elem) == CNhbr.second.end()){
            stdev += pow((norm(cc.meanNhbr[CNhbr.first] - elem->pos)),2);
            countele++; 
          }
        }
        forall(const vertex &elem, T.neighbors(CNhbr.first.second)){
          if(CNhbr.second.find(elem) == CNhbr.second.end()){
            stdev += pow((norm(cc.meanNhbr[CNhbr.first] - elem->pos)),2);
            countele++; 
          }
        }
        stdev /= countele;
        stdev = pow(stdev, 0.5); 
        if(stdev < cc.StDev[CNhbr.first.first->label])
          cc.StDev[CNhbr.first.first->label] = stdev;
        if(stdev < cc.StDev[CNhbr.first.second->label])
          cc.StDev[CNhbr.first.second->label] = stdev;
    }

    //Compute closest neighbour map
    forall(const vertex &v, T){
      forall(const vertex &n, cc.Nhbr[v->label]){  
        if(cc.computed.find(v) == cc.computed.end() and cc.StDev[v->label] == cc.StDev[n->label] and cc.NhbrElem[std::make_pair(v, n)].size() > 1 ){
          cc.computed.insert(v);
          cc.computed.insert(n);
          cc.NhbrPr[v] = n;
          cc.NhbrPr[n] = v;
        } 
      }
    }
    forall(VertexPr Pr, cc.NhbrPr){
      cc.NeighbourMap.insert(std::make_pair(Pr.first, Pr.second));
      std::cout<<"\n"<<Pr.first->label<<"\t"<<Pr.second->label;
    }
    std::cout<<"\n cc.NeighbourMap size\t"<< cc.NeighbourMap.size();
    cc.FeatureVector.push_back(cc.StDev);
  }
  void CalculateCommonNhbrs(vvGraph &T, CellCluster &cc){
    forall(VertexPr Pr, cc.NeighbourMap){
      int RepeatedNeighbours = 0;
      forall(const vertex &elem, cc.Nhbr[Pr.first->label]){
        if(cc.Nhbr[Pr.second->label].find(elem) != cc.Nhbr[Pr.second->label].end()){
          RepeatedNeighbours++; 
        }
      }
      cc.TotalNhbr[Pr.first->label] = cc.TotalNeighbours[Pr.first->label] + cc.TotalNeighbours[Pr.second->label] - RepeatedNeighbours - 2;
      cc.TotalNhbr[Pr.second->label] = cc.TotalNeighbours[Pr.first->label] + cc.TotalNeighbours[Pr.second->label] - RepeatedNeighbours- 2;
      std::cout<<"\nTotal nhbrs\t"<<cc.TotalNeighbours[Pr.first->label]<<"\t"<<cc.TotalNeighbours[Pr.second->label]<<"\t"<<RepeatedNeighbours<<"\t"<<cc.TotalNhbr[Pr.first->label];
    } 
    cc.FeatureVector.push_back(cc.TotalNhbr);
  }
  void CalculateAverage(vvGraph& T, CellCluster &cc){
    std::cout<<"\n\t"<<"v->label"<<"\t"<<"cc.Area[v->label]"<<"\t"<<"cc.Perimeter[v->label]"<<"\t"
      <<"cc.stomataBending[v->label]"<<"\t"<<"cc.DiffLenBr[v->label]"<<"\t"<<"cc.Bending[v->label]"<<"\t"<<"cc.TotalNhbr[v->label]"<<"\t"
      <<"cc.TotalNeighbours[v->label]"<<"\t"<<"cc.StDev[v->label]"<<"\t"<<"\t"<<"cc.Growth[v->label]"<<"\t"
      <<"cc.StrMax[v->label]"<<"\t"<<"cc.StrMin[v->label]"<<"\t"<<"cc.Aniso[v->label]"<<"\t"<<"cc.Curvature[v->label]";
    forall(const vertex v, T){
      if(v->label >0)
      std::cout<<"\n\t"<<v->label<<"\t"<<cc.Area[v->label]<<"\t"<<cc.Perimeter[v->label]<<"\t"
      <<cc.stomataBending[v->label]<<"\t"<<cc.DiffLenBr[v->label]<<"\t"<<cc.Bending[v->label]<<"\t"<<cc.TotalNhbr[v->label]<<"\t"
      <<cc.TotalNeighbours[v->label]<<"\t"<<cc.StDev[v->label]<<"\t"<<"\t"<<cc.Curvature[v->label];
    }
    forall(const vertex &v, T){
      if(v->label != -1){
        cc.AvgArea[v->label] = 0;
        cc.AvgPerimeter[v->label] = 0;
        cc.AvgstomataBending[v->label] = 0;
        cc.AvgDiffLenBr[v->label] = 0;
        cc.AvgBending[v->label]  = 0;
        cc.AvgTotalNhbr[v->label] = 0;
        cc.AvgStandardDeviation[v->label] = 0;
        cc.AvgTotalNeighbours[v->label]  = 0;
        cc.AvgCurvature[v->label]  = 0;
        forall(const vertex &n, cc.Nhbr[v->label]){
          cc.AvgArea[v->label] += cc.Area[n->label];
          cc.AvgPerimeter[v->label] += cc.Perimeter[n->label];
          cc.AvgstomataBending[v->label] += cc.stomataBending[n->label];
          cc.AvgDiffLenBr[v->label] += cc.DiffLenBr[n->label];
          cc.AvgBending[v->label] += cc.Bending[n->label];
          cc.AvgTotalNhbr[v->label] += cc.TotalNhbr[n->label];
          cc.AvgStandardDeviation[v->label] += cc.StDev[n->label];
          cc.AvgTotalNeighbours[v->label] += cc.TotalNeighbours[n->label];
          cc.AvgCurvature[v->label] += cc.Curvature[n->label];
        }
        cc.AvgArea[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgPerimeter[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgstomataBending[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgDiffLenBr[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgBending[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgTotalNhbr[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgStandardDeviation[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgTotalNeighbours[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgCurvature[v->label] /= cc.Nhbr[v->label].size();
      }
    }
    cc.FeatureVector.push_back(cc.AvgArea);
    cc.FeatureVector.push_back(cc.AvgPerimeter);
    cc.FeatureVector.push_back(cc.AvgstomataBending);
    cc.FeatureVector.push_back(cc.AvgDiffLenBr);
    cc.FeatureVector.push_back(cc.AvgBending);
    cc.FeatureVector.push_back(cc.AvgTotalNhbr);
    cc.FeatureVector.push_back(cc.AvgStandardDeviation);
    cc.FeatureVector.push_back(cc.AvgTotalNeighbours);
    cc.FeatureVector.push_back(cc.AvgCurvature);    
  }
  void FeatureScaling(std::vector<borderType> FeatureVector){ 
    forall(borderType feature, FeatureVector){
      float min = HUGE_VAL, max = 0, sum = 0;
      forall(to_loop_map loop1, feature){
        if(max < loop1.second)
          max = loop1.second;
        if(min > loop1.second)
          min = loop1.second;
        sum += loop1.second;
      }
      sum /= feature.size();
      forall(to_loop_map loop1, feature){
        std::cout<<"b "<<loop1.second;
        loop1.second = (sum - loop1.second)/(max - min);
        if(loop1.second < 0)
          loop1.second *= -1;
        std::cout<<"a "<<loop1.second<<"\n";
      }
    }
  }
  void ComputeData(vvGraph &T, CellCluster &cc){
    CalculateArea(T, cc);
    CalculatePerimeter(T, cc);
    CalculateBending(T, cc);
    CalculateLenBr(T, cc);
    CalculateNeighbours(T, cc);
    CalculateCurvature(T, cc);
    FindCommonNeighbours(T, cc);
    CalculateStandardDeviation(T, cc);
    CalculateStomatalBending(T, cc); 
    CalculateCommonNhbrs(T, cc);
    CalculateAverage(T, cc);
    FeatureScaling(cc.FeatureVector);
    cc.FeatureMap[cc.Area] = "Area";
    cc.FeatureMap[cc.Perimeter] = "Perimeter";
    cc.FeatureMap[cc.Bending] = "Bending";
    cc.FeatureMap[cc.stomataBending] = "Stomata Bending";
    cc.FeatureMap[cc.DiffLenBr] = "Length Breadth Ratio";
    cc.FeatureMap[cc.TotalNhbr] = "Total Nhbrs of 2 similar cells";
    cc.FeatureMap[cc.TotalNeighbours] = "Total Neighbours";
    cc.FeatureMap[cc.Curvature] = "Curvature";
    cc.FeatureMap[cc.StDev] = "Standard Deviation";
    cc.FeatureMap[cc.AvgArea] = "Avg. Area";
    cc.FeatureMap[cc.AvgPerimeter] = "Avg. Perimeter";
    cc.FeatureMap[cc.AvgstomataBending] = "Avg. Stomatal Bending";
    cc.FeatureMap[cc.AvgBending] = "Avg. Bending";
    cc.FeatureMap[cc.AvgTotalNhbr] = "Avg. Total Neighbours of similar cells";
    cc.FeatureMap[cc.AvgStandardDeviation] = "Avg. Standard Deviation";
    cc.FeatureMap[cc.AvgTotalNeighbours] = "Avg. Total Neighbours";
    cc.FeatureMap[cc.AvgCurvature] = "Avg. Curvature";
    
  }
  void LoadValues(vvGraph& T,QString filename1, QString filename2, QString filename3, QString filename4,CellCluster &cc){
    std::map<int, double> GrowthTemp, StrMaxTemp, StrMinTemp, AnisoTemp;
    typedef std::pair<int,double> LoopGrwth;
    std::cout<<"\n"<<filename1.toLocal8Bit().constData()<<"\t"<<filename2.toLocal8Bit().constData()<<"\t"<<filename3.toLocal8Bit().constData()<<"\t"<<filename4.toLocal8Bit().constData();
    if(filename1 != "Growth Data"){
      QFile file(filename1);
      if(!file.open(QIODevice::ReadOnly))
        {
          std::cout<<"\nWrong";
//        setErrorMessage(QString("File '%1' cannot be opened for reading").arg(filename));
          //return false;
      }
      QTextStream ss(&file);
      QString line = ss.readLine();
      QStringList fields = line.split(",");
      // variables to be filled
      bool ok;
      std::cout << "start" << ss.status();
      while(ss.status() == QTextStream::Ok){
        std::cout<<"\nIn while True";
        fields = ss.readLine().split(",");
        std::cout << "\nstatus\t" << ss.status();
        if(fields[0].isEmpty()) break;
        int label = fields[0].toInt(&ok);
        double value = fields[1].toDouble(&ok);
        std::cout << "\nlll  " << label << " value " << value <<"\t feild 1"<<fields[0].toDouble(&ok);
        //cc.Growth[label] = value;
        GrowthTemp[label] = value;
      }
    }
    forall(LoopGrwth loop, GrowthTemp)
      std::cout<<"\nMap1\t"<<loop.first<<"\t"<<loop.second;
    if(filename2 != "Stretch Max Data"){
      QFile file1(filename2);
      if(!file1.open(QIODevice::ReadOnly))
        {
          std::cout<<"\nWrong";
      }
      QTextStream ss1(&file1);
      QString line = ss1.readLine();
      QStringList fields = line.split(",");
      // variables to be filled
      bool ok;
      //cout << "start" << endl;
      while(ss1.status() == QTextStream::Ok){
        fields = ss1.readLine().split(",");
        if(fields[0].isEmpty()) break;
        int label = fields[0].toInt(&ok);
        double value = fields[1].toDouble(&ok);
        std::cout << "\nlll  " << label << " value " << value;
        //cc.StrMax[label] = value;
        StrMaxTemp[label] = value;
      }
    }
    forall(LoopGrwth loop, StrMaxTemp)
      std::cout<<"\nMap2\t"<<loop.first<<"\t"<<loop.second;
    if(filename3 != "Stretch Min Data"){
      QFile file2(filename3);
      if(!file2.open(QIODevice::ReadOnly))
        {
          std::cout<<"\nWrong";
      }
      QTextStream ss2(&file2);
      QString line = ss2.readLine();
      QStringList fields = line.split(",");
      // variables to be filled
      bool ok;
      //cout << "start" << endl;
      while(ss2.status() == QTextStream::Ok){
        fields = ss2.readLine().split(",");
        if(fields[0].isEmpty()) break;
        int label = fields[0].toInt(&ok);
        double value = fields[1].toDouble(&ok);
        std::cout << "\nlll  " << label << " value " << value;
        //cc.StrMin[label] = value;
        StrMinTemp[label] = value;
      }
    }
    forall(LoopGrwth loop, StrMinTemp)
      std::cout<<"\nMap3\t"<<loop.first<<"\t"<<loop.second;
    if(filename4 != "Anisotropy Data"){
      QFile file3(filename4);
      if(!file3.open(QIODevice::ReadOnly))
        {
          std::cout<<"\nWrong";
      }
      QTextStream ss3(&file3);
      QString line = ss3.readLine();
      QStringList fields = line.split(",");
      // variables to be filled
      bool ok;
      //cout << "start" << endl;
      while(ss3.status() == QTextStream::Ok){
        fields = ss3.readLine().split(",");
        if(fields[0].isEmpty()) break;
        int label = fields[0].toInt(&ok);
        double value = fields[1].toDouble(&ok);
        std::cout << "\nlll  " << label << " value " << value;
        //cc.Aniso[label] = value;
        AnisoTemp[label] = value;
      }
    }
    forall(LoopGrwth loop, AnisoTemp)
      std::cout<<"\nMap4\t"<<loop.first<<"\t"<<loop.second;

    cc.Growth = GrowthTemp;
    cc.StrMax = StrMaxTemp;
    cc.StrMin = StrMinTemp;
    cc.Aniso = AnisoTemp;
    forall(const vertex &v, T){
      if(v->label != -1){
        cc.AvgGrowth[v->label]  = 0;
        cc.AvgStrMax[v->label]  = 0;
        cc.AvgStrMin[v->label]  = 0;
        cc.AvgAniso[v->label]  = 0;
        forall(const vertex &n, cc.Nhbr[v->label]){
          cc.AvgGrowth[v->label] += cc.Growth[n->label];
          cc.AvgStrMax[v->label] += cc.StrMax[n->label];
          cc.AvgStrMin[v->label] += cc.StrMin[n->label];
          cc.AvgAniso[v->label] += cc.Aniso[n->label];
        }
        cc.AvgGrowth[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgStrMax[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgStrMin[v->label] /= cc.Nhbr[v->label].size();
        cc.AvgAniso[v->label] /= cc.Nhbr[v->label].size();
      }
    }
    cc.FeatureMap[cc.Growth] = "Growth";
    cc.FeatureMap[cc.StrMax] = "Str Maximum";
    cc.FeatureMap[cc.StrMin] = "Str Minimum";
    cc.FeatureMap[cc.Aniso] = "Anisotropy";
    cc.FeatureMap[cc.AvgGrowth] = "Avg. Growth";
    cc.FeatureMap[cc.AvgStrMax] = "Avg. StrMax";
    cc.FeatureMap[cc.AvgStrMin] = "Avg. StrMin";
    cc.FeatureMap[cc.AvgAniso] = "Avg. Anisotropy";
  }     

  bool SizeCluster_SVM_Prepare::run(Mesh *mesh, QString filename1, QString filename2, QString filename3, QString filename4, CellCluster &cc)
  {  
    // Fix Corners
    FixMeshCorners FixCorners(*this);
    FixCorners.run(mesh, false, true, 0);
    std::cout<<__LINE__<<"\n";
    // Make cells
    bool result; 
    vvGraph T;
    mgxmToMgxc(mesh->graph(), T, 0);
    ConvertToMgxc cnvrt(*this);
    cnvrt.run(mesh, 1);
    if(T.size() == 0)  //If cell mesh is present 
      T = mesh->graph();
    ComputeData(T, cc);  //Compute data
    std::cout<<__LINE__<<"\n";
    //LoadValues(T, filename1, filename2, filename3, filename4, cc);
    return true;
  }
  REGISTER_PROCESS(SizeCluster_SVM_Prepare);

  bool SizeCluster_SVM_DisplayCellData::FeatureBox(CellCluster &cc){  
    forall(TypeFeatureMap Feature, cc.FeatureMap) {
      QTreeWidgetItem *item = new QTreeWidgetItem(QStringList() << Feature.second);
      item->setFlags(Qt::ItemIsEnabled | Qt::ItemIsSelectable | Qt::ItemIsUserCheckable);
      item->setCheckState(1, Qt::Unchecked);
      ui.svmTreeWidget->addTopLevelItem(item);
    }
    ui.svmTreeWidget->sortItems(0, Qt::AscendingOrder);
    for(int i = 0; i < 2; i++)
      ui.svmTreeWidget->resizeColumnToContents(i);
  }

  bool SizeCluster_SVM_DisplayCellData::initialize(QStringList &parms, QWidget *parent){ 
    Mesh *m1 = currentMesh();
    m1->showHeat();
    vvGraph &T = m1->graph(); 
    std::cout<<"\nGraph size\t"<<T.size()<<"\n";
    QDialog dlg(parent);
    ui.setupUi(&dlg);
    FeatureBox(cc);
    connect(ui.svmTreeWidget, SIGNAL(itemClicked(QTreeWidgetItem *, int)),
            this, SLOT(on_svmTreeWidget_itemClicked(QTreeWidgetItem *, int )));
    connect(ui.selectAllButton, SIGNAL(clicked()), this, SLOT(on_selectAllButton_clicked()));
    dlg.exec();
    selectFeatures();
    std::cout<<__LINE__ <<"\n";
  }
  void SizeCluster_SVM_DisplayCellData::on_svmTreeWidget_itemClicked(QTreeWidgetItem *item, int column){
    Mesh *m1 = currentMesh();
    m1->showHeat();
    std::map <int, double> heatMap;
    float maxValue = 0, minValue = HUGE_VAL;
    vvGraph &T = m1->graph(); 
    std::cout<<"\nHeat map\t"<< heatMap.size() << "\t" << cc.FeatureMap.size()<<"\t"<<T.size();
    forall(TypeFeatureMap Feature, cc.FeatureMap){
      std::cout<<__LINE__ <<"\n";
      if(item->text(0) == Feature.second){
        std::cout<<__LINE__ <<"\n";
        heatMap = Feature.first;
        qDebug() << item->text(0);
        std::cout<<"\nHeat map in loop\t"<< heatMap.size() <<"\t"<< Feature.first.size()<<"\n";
        break;
      }
      std::cout<<__LINE__ <<"\n";
    }
    std::cout<<"\nAfter Heat map\t"<< heatMap.size()<<"\n";
    if(heatMap.size() > 0){
      std::cout<<__LINE__ <<"\n";
      forall(to_loop_map loop, heatMap){
        if(minValue > heatMap[loop.first])
          minValue = heatMap[loop.first];
        if(maxValue < heatMap[loop.first])
          maxValue = heatMap[loop.first];
      }
      std::cout<<__LINE__ <<"\n";
      m1->wallHeat().clear();
      m1->labelHeat().clear();
      forall(const vertex v, T){
        if(v->label != -1)
          m1->labelHeat()[v->label] = heatMap[v->label];
      }
      std::cout<<__LINE__ <<"\n";
      m1->heatMapBounds() = Point2f(minValue, maxValue);
      m1->updateTriangles();
      m1->showHeat();
      
      //m1->setShowLabel("Label Heat");
      m1->updateAll();
      std::cout<<__LINE__ <<"\n";
      updateState();
      updateViewer();
    }
  }
  void SizeCluster_SVM_DisplayCellData::on_selectAllButton_clicked(){
    for(int itemCount = 0; itemCount < ui.svmTreeWidget->topLevelItemCount(); ++itemCount){
      ui.svmTreeWidget->topLevelItem(itemCount)->setCheckState(1, Qt::Checked);
      qDebug()<<"\t"<<ui.svmTreeWidget->topLevelItem(itemCount)->text(0)<<"\t"<<ui.svmTreeWidget->topLevelItem(itemCount)->checkState(1);
    }
  }
  void SizeCluster_SVM_DisplayCellData::selectFeatures(){ 
    std::vector<QString> Features;
    cc.SelectedFeatures.clear();
    for(int itemCount = 0; itemCount < ui.svmTreeWidget->topLevelItemCount(); ++itemCount){
      std::cout<<__LINE__<<"\n"; 
      qDebug()<<"\t"<<ui.svmTreeWidget->topLevelItem(itemCount)->checkState(1);
      if(Qt::Checked == ui.svmTreeWidget->topLevelItem(itemCount)->checkState(1)){
          std::cout<<__LINE__<<"\n";
          //qDebug()<<"\t"<<ui.svmTreeWidget->topLevelItem(itemCount)->text(0);
          Features.push_back(ui.svmTreeWidget->topLevelItem(itemCount)->text(0));
        }
    }
    std::cout<<"\n"<<Features.size()<<"\t"<<__LINE__<<"\n";
    forall(TypeFeatureMap Feature, cc.FeatureMap){
      forall(QString f, Features){
        //qDebug()<<"\nFeature\t"<<f;
        if(Feature.second == f)
          cc.SelectedFeatures.push_back(Feature.first);
      }
    }
    std::cout<<__LINE__<<"\n";
    std::cout<<"Selected feaures\t"<<cc.SelectedFeatures.size()<<"\n";
  }
  REGISTER_PROCESS(SizeCluster_SVM_DisplayCellData);

  bool SizeCluster_SVM::run(Mesh *mesh, QString choice, CellCluster &cc)
   {  
      vvGraph& S = mesh->graph();
      if(cc.flag == 1){
        forall(const vertex& v, S){
          cc.countlabel[v->label] = 0; 
        }
      }
      cc.flag ++;
      int i = 0;
      //For cohort of selected cells 
      if (choice == "Giant Cells"){
        forall(const vertex &v, S){
          if(cc.countlabel[v->label]<1 && v->selected && v->label>0){
            sample_type samp;
            i = 0;
            forall(borderType map, cc.SelectedFeatures){
              //std::cout<<"\nfor i \t"<<i;
              forall(to_loop_map loop,map){
                label = loop.first;
                if(label == v->label){
                  samp(i) = loop.second;
                  //std::cout<<"\t"<<i<<"\t"<<loop.first<<"\t"<<samp(i);
                  i++;
                }
              }
            }
            cc.samples.push_back(samp);
            cc.labels.push_back(1);
            mesh->parents()[v->label] = 1;
            cc.countlabel[v->label]++;
          }
        }
        cc.x++;
      }
      else if (choice == "Old Stomata"){
        forall(const vertex &v, S){
          if(cc.countlabel[v->label]<1 && v->selected && v->label>0){
            sample_type samp;
            i = 0;
            forall(borderType map, cc.SelectedFeatures){
              forall(to_loop_map loop,map){
                label = loop.first;
                if(label == v->label){
                  samp(i) = loop.second;
                  //std::cout<<"\t"<<i<<"\t"<<loop.first<<"\t"<<samp(i);
                  i++;
                }
              }
            }
            cc.samples.push_back(samp);
            cc.labels.push_back(2);
            mesh->parents()[v->label] = 2;
            cc.countlabel[v->label]++;
          }
        }
        cc.x++;
      }
      else if (choice == "Young Stomata"){
        forall(const vertex &v, S){
          if(cc.countlabel[v->label] < 1 && v->selected && v->label > 0){
            sample_type samp;
            i = 0;
            forall(borderType map, cc.SelectedFeatures){
              std::cout<<"\nfor i \t"<<i;
              forall(to_loop_map loop,map){
                label = loop.first;
                if(label == v->label){
                  samp(i) = loop.second;
                  std::cout<<"\t"<<i<<"\t"<<loop.first<<"\t"<<samp(i);
                  i++;
                }
              }
            }
            cc.samples.push_back(samp);
            cc.labels.push_back(3);
            mesh->parents()[v->label] = 3;
            cc.countlabel[v->label]++;
          }
        }
        cc.x++;
      }
      else if (choice == "Large Regular Cells"){
        forall(const vertex &v, S){
          if(cc.countlabel[v->label]<1 && v->selected && v->label>0){
            sample_type samp;
            i = 0;
            forall(borderType map,cc.SelectedFeatures){
              std::cout<<"\nfor i \t"<<i;
              forall(to_loop_map loop,map){
                label = loop.first;
                if(label == v->label){
                  samp(i) = loop.second;
                  std::cout<<"\t"<<i<<"\t"<<loop.first<<"\t"<<samp(i);
                  i++;
                }
              }
            }
            cc.samples.push_back(samp);
            cc.labels.push_back(4);
            mesh->parents()[v->label] = 4;
            cc.countlabel[v->label]++;
          }
        }
        cc.x++;
      }
      else if (choice == "Small Regular Cells"){
        forall(const vertex &v, S){
          if(cc.countlabel[v->label]<1 && v->selected && v->label>0){
            sample_type samp;
            i = 0;
            forall(borderType map,cc.SelectedFeatures){
              std::cout<<"\nfor i \t"<<i;
              forall(to_loop_map loop,map){
                label = loop.first;
                if(label == v->label){
                  samp(i) = loop.second;
                  std::cout<<"\t"<<i<<"\t"<<loop.first<<"\t"<<samp(i);
                  i++;
                }
              }
            }
            cc.samples.push_back(samp);
            cc.labels.push_back(5);
            mesh->parents()[v->label] = 5;
            cc.countlabel[v->label]++;
          }
        }
        cc.x++;
      }
      else if (choice == "Large Marginal Cells"){
        forall(const vertex &v, S){
          if(cc.countlabel[v->label]<1 && v->selected && v->label>0){
            sample_type samp;
            i = 0;
            forall(borderType map,cc.SelectedFeatures){
              std::cout<<"\nfor i \t"<<i;
              forall(to_loop_map loop,map){
                label = loop.first;
                if(label == v->label){
                  samp(i) = loop.second;
                  std::cout<<"\t"<<i<<"\t"<<loop.first<<"\t"<<samp(i);
                  i++;
                }
              }
            }
            cc.samples.push_back(samp);
            cc.labels.push_back(6);
            mesh->parents()[v->label] = 6;
            cc.countlabel[v->label]++;
          }
        }
        cc.x++;
      }
      else if (choice == "Small Marginal Cells"){
        forall(const vertex &v, S){
          if(cc.countlabel[v->label]<1 && v->selected && v->label>0){
            sample_type samp;
            i = 0;
            forall(borderType map,cc.SelectedFeatures){
              std::cout<<"\nfor i \t"<<i;
              forall(to_loop_map loop,map){
                label = loop.first;
                if(label == v->label){
                  samp(i) = loop.second;
                  std::cout<<"\t"<<i<<"\t"<<loop.first<<"\t"<<samp(i);
                  i++;
                }
              }
            }
            cc.samples.push_back(samp);
            cc.labels.push_back(7);
            mesh->parents()[v->label] = 7;    
            cc.countlabel[v->label]++;
          }
        }
        cc.x++;
      }
      else if (choice == "Petiole Cells"){
        forall(const vertex &v, S){
          if(cc.countlabel[v->label]<1 && v->selected && v->label>0){
            sample_type samp;
            i = 0;
            forall(borderType map,cc.SelectedFeatures){
              std::cout<<"\nfor i \t"<<i;
              forall(to_loop_map loop,map){
                label = loop.first;
                if(label == v->label){
                  samp(i) = loop.second;
                  std::cout<<"\t"<<i<<"\t"<<loop.first<<"\t"<<samp(i);
                  i++;
                }
              }
            }
            cc.samples.push_back(samp);
            cc.labels.push_back(8);
            mesh->parents()[v->label] = 8;
            cc.countlabel[v->label]++;
          }
        }
        cc.x++;
      }
      else if (choice == "Pavement Cells"){
        forall(const vertex &v, S){
          if(cc.countlabel[v->label]<1 && v->selected && v->label>0){
            sample_type samp;
            i = 0;
            forall(borderType map,cc.SelectedFeatures){
              std::cout<<"\nfor i \t"<<i;
              forall(to_loop_map loop,map){
                label = loop.first;
                if(label == v->label){
                  samp(i) = loop.second;
                  std::cout<<"\t"<<i<<"\t"<<loop.first<<"\t"<<samp(i);
                  i++;
                }
              }
            }
            cc.samples.push_back(samp);
            cc.labels.push_back(9);
            mesh->parents()[v->label] = 9;
            cc.countlabel[v->label]++;
          }
        }
        cc.x++;
      }
      mesh->updateTriangles();
    return true;
  }
  REGISTER_PROCESS(SizeCluster_SVM);
  bool SizeCluster_SVM_Clear::run(Mesh *mesh, CellCluster &cc){  
      CellCluster ccNew;
      vec labels;
      std::vector<sample_type> samples;
      std::vector<borderType> M;
      cc = ccNew;
      cc.SelectedFeatures = M;
      cc.samples = samples;
      cc.labels = labels;
      cc.x = 0;
      vvGraph& S = mesh->graph();
      forall(const vertex& v, S){
        cc.countlabel[v->label] = 0; 
        mesh->parents()[v->label] = 0;
      }

      mesh->updateAll();
      return true;
   }
   REGISTER_PROCESS(SizeCluster_SVM_Clear);

   bool SizeCluster_SVM_Train::run(Mesh *mesh, QString choice, CellCluster &cc){ 
      vvGraph& S = mesh->graph();
      if(choice == "Train using SVM"){

        }
      }  
      if(choice == "Train using Perceptron"){
         
      }
      return true;
    }
  REGISTER_PROCESS(SizeCluster_SVM_Train);
         
  bool SizeCluster_SVM_Test::run(Mesh *mesh, QString choice, CellCluster &cc)
    { 
      deserialize("function.dat") >> df3;
      vvGraph& S = mesh->graph();
      std::map<int,int> count;
      typedef std::pair<int,double> to_loop_borderType;
      if (choice == "Classify Cells"){
        forall(const vertex &v, S){
      }
      if (choice == "All cells"){
        forall(pr Pr, cc.LabelParentMap){
          mesh->parents()[Pr.first] = Pr.second;
        }
      }    
      mesh->updateTriangles();
      return true;
    }
   REGISTER_PROCESS(SizeCluster_SVM_Test);

   bool SizeCluster_SVM_Post::run(Mesh *mesh, QString ch1, QString ch2)
    { 
      std::map<int, int> MapChange;
      typedef std::pair<int, int> pr;
      vvGraph& S = mesh->graph();
      if(ch1 > 0){
        forall(const vertex &v, S){
          if(v->selected)
            mesh->parents()[v->label] = ch1.toInt();  
        }
      }
      if(ch2 == "Yes"){
        forall(const vertex &v, S){
          if(mesh->parents()[v->label] == 4)
            mesh->parents()[v->label] = 5;  
        } 
      }

      mesh->updateAll();
      return true;
    }
   REGISTER_PROCESS(SizeCluster_SVM_Post);

   bool SizeCluster_SVM_Calc::run(QString ch1, QString ch2)
    { 
    std::map<int, int> Map1, Map2;
    typedef std::pair<int, int> pr;
    float cor1 = 0, cor2 = 0,cor3 = 0, cor4 = 0,cor5 = 0, cor6 = 0,cor7 = 0, cor8 = 0, cor9 = 0, totcor = 0, tot = 0;
    
    if(ch1 != " " and ch2 != " "){
      QFile file(ch1);
      if(!file.open(QIODevice::ReadOnly))
        {
          std::cout<<"\nWrong";
//        setErrorMessage(QString("File '%1' cannot be opened for reading").arg(filename));
          //return false;
      }
      QTextStream ss(&file);
      QString line = ss.readLine();
      QStringList fields = line.split(",");
      // variables to be filled
      bool ok;
      std::cout << "start" << ss.status();
      while(ss.status() == QTextStream::Ok){
        std::cout<<"\nIn while True";
        fields = ss.readLine().split(",");
        std::cout << "\nstatus\t" << ss.status();
        if(fields[0].isEmpty()) break;
        int label1 = fields[0].toInt(&ok);
        int parLabel1 = fields[1].toInt(&ok);
        std::cout << "\nlabel " << label1 << " Parent label " << parLabel1;
        Map1[label1] = parLabel1;
      }
      QFile file1(ch2);
      if(!file1.open(QIODevice::ReadOnly))
        {
          std::cout<<"\nWrong";
      }
      QTextStream ss1(&file1);
      line = ss1.readLine();
      fields = line.split(",");
      while(ss1.status() == QTextStream::Ok){
        fields = ss1.readLine().split(",");
        if(fields[0].isEmpty()) break;
        int label2 = fields[0].toInt(&ok);
        int parLabel2 = fields[2].toInt(&ok);
        std::cout << "\nlabel " << label2 << " Parent label " << parLabel2;
        Map2[label2] = parLabel2;
      }
    }  
    forall(pr Pr1, Map1){
      forall(pr Pr2, Map2){
        if(Pr1.first == Pr2.first){
          if(Pr1.second == Pr2.second and Pr1.second == 1)
            cor1 ++;
          if(Pr1.second == Pr2.second and Pr1.second == 2)
            cor2 ++;
          if(Pr1.second == Pr2.second and Pr1.second == 3)
            cor3 ++;
          if(Pr1.second == Pr2.second and Pr1.second == 4)
            cor4 ++;
          if(Pr1.second == Pr2.second and Pr1.second == 5)
            cor5 ++;
          if(Pr1.second == Pr2.second and Pr1.second == 6)
            cor6 ++;
          if(Pr1.second == Pr2.second and Pr1.second == 7)
            cor7 ++;
          if(Pr1.second == Pr2.second and Pr1.second == 8)
            cor8 ++;
          if(Pr1.second == Pr2.second and Pr1.second == 9)
            cor9 ++;
          tot++;
        }
      }
    }
    totcor = (cor1 + cor2 + cor3 + cor4 + cor5 + cor6 + cor7 + cor8 + cor9)/tot;
    std::cout<<"\nAccuracy on entire mesh\t"<<totcor;  
    std::cout<<"\nAccuracy on cell type 1 - Giant cells\t"<<(cor1/tot);    
    std::cout<<"\nAccuracy on cell type 2 - Old Stomata\t"<<(cor2/tot); 
    std::cout<<"\nAccuracy on cell type 3 - Young Stomata\t"<<(cor3/tot); 
    std::cout<<"\nAccuracy on cell type 4 - Large Regular Cells\t"<<(cor4/tot); 
    std::cout<<"\nAccuracy on cell type 5 - Small Regular Cells\t"<<(cor5/tot); 
    std::cout<<"\nAccuracy on cell type 6 - Large Marginal Cells\t"<<(cor6/tot); 
    std::cout<<"\nAccuracy on cell type 7 - Small Marginal Cells\t"<<(cor7/tot); 
    std::cout<<"\nAccuracy on cell type 8 - Petiole Cells\t"<<(cor8/tot); 
    std::cout<<"\nAccuracy on cell type 9 - Pavement Cells\t"<<(cor9/tot); 
      return true;
    
    }
   REGISTER_PROCESS(SizeCluster_SVM_Calc);



 // file dialogue for save data
 bool SaveSVMData::initialize(QStringList& parms, QWidget* parent)
  {
    QString &filename = parms[0];
    if(filename.isEmpty())
      filename = QFileDialog::getSaveFileName(0, "Choose spreadsheet file to save", QDir::currentPath(), "CSV files (*.csv)");
    if(filename.isEmpty())
      return false;
    if(!filename.endsWith(".csv", Qt::CaseInsensitive))
      filename += ".csv";
    parms[0] = filename;
    return true;
  }


  bool SaveSVMData::run(QString filename, Mesh *m1, CellCluster &cc)
  {
    //cout << "file   " << filename.toLocal8Bit().constData() << endl;

    // file for root cell data
    QFile file(filename);
    if(!file.open(QIODevice::WriteOnly))
    {
      setErrorMessage(QString("File '%1' cannot be opened for writing").arg(filename));
      return false;
    }
    QTextStream out(&file);

    out << "Label,Parent,Cell Area,Perimeter,Stomatal Bending,Len. Br Ratio,Bending,Total neighbours incl. closest neighbour,Avg Radius,Total Neighbours single cell,Growth,Stretch Max,Stretch Min,Anisotropy,PeriSq/Area,Avg Cell Area,Avg Perimeter,Avg Stomatal Bending,Avg Len. Br Ratio,Avg Bending,Avg Total neighbours incl. closest neighbour,Avg Avg Radius,Avg Total Neighbours single cell,Avg Growth,Avg Stretch Max,Avg Stretch Min,Avg Anisotropy,Avg PeriSq/Area" << endl;

   // Mesh *m1 = mesh(0);
    vvGraph& S = m1->graph();

    forall(const pr &p, cc.Area){
      int l = p.first;

    out << l << "," << m1->parents()[l] << ","  << cc.Area[l] << "," << cc.Perimeter[l] << "," << cc.stomataBending[l] << "," << cc.DiffLenBr[l] << "," << cc.Bending[l] << "," << cc.TotalNhbr[l] << "," << cc.AvgRadius[l] << "," << cc.TotalNeighbours[l] << "," << cc.Growth[l] << "," << cc.StrMax[l] << "," << cc.StrMin[l] << "," << cc.Aniso[l] << "," << cc.Curvature[l] << "," << cc.AvgArea[l] << "," << cc.AvgPerimeter[l] << "," << cc.AvgstomataBending[l] << "," << cc.AvgDiffLenBr[l] << "," << cc.AvgBending[l] << "," << cc.AvgTotalNhbr[l] << "," << cc.AvgAvgRadius[l] << "," << cc.AvgTotalNeighbours[l] << "," << cc.AvgGrowth[l] << "," << cc.AvgStrMax[l] << "," << cc.AvgStrMin[l] << "," << cc.AvgAniso[l] << "," << cc.AvgCurvature[l] << endl;

    }

    return true;
  }
  REGISTER_PROCESS(SaveSVMData);

}

 
