//
// This file is part of MorphoGraphX - http://www.MorphoGraphX.org
// Copyright (C) 2012-2015 Richard S. Smith and collaborators.
//
// If you use MorphoGraphX in your work, please cite:
//   http://dx.doi.org/10.7554/eLife.05864
//
// MorphoGraphX is free software, and is licensed under under the terms of the 
// GNU General (GPL) Public License version 2.0, http://www.gnu.org/licenses.
//
#include <SizeCluster_SVM.hpp>
#include <CellCluster.hpp>
#include <math.h> 
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_vector.h>
#include <gsl/gsl_eigen.h>
#include <iostream>
#include <vector>
#include <fstream>
#include <limits>
#include <Progress.hpp>
#include <QFileDialog>
#include <QFileDialog>


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
    #define Malloc(type,n) (type *)malloc((n)*sizeof(type))
    int nr_fold;
    
    bool double_equals(double a, double b, double epsilon = 0.001){
        return std::abs(a - b) < epsilon;
    }

    bool SizeCluster_SVM_Prepare::FeatureScaling(vvGraph &T, IntFloatAttr &map){ 
        float min = HUGE_VAL, max = 0, sum = 0, standDev = 0;
        forall(to_loop_map loop1, map){
            if(loop1.first == -1)
                continue;
            //std::cout<<loop1.first<<" "<<loop1.second<<"\n";
            sum += loop1.second;
            if(max < loop1.second)
              max = loop1.second;
            if(min > loop1.second)
              min = loop1.second;
        }
        std::cout<<"\nmax\t"<<max<<"\tmin\t"<<min<<"\tsum\t"<<sum<<"\t"<<map.size();
        sum /= map.size();
        //TODO: Find std dev and fix scaling.
        //Did : Changed FeatureVector, FeatureMap in terms of keys and subsequently the initialization window.
        //Calculate Std dev
        forall(to_loop_map loop1, map){
          standDev += pow((loop1.second - sum), 2);
        }
        standDev /= map.size();
        standDev = pow(standDev, 0.5); 
        std::cout<<"\nmax\t"<<max<<"\tmin\t"<<min<<"\tsum\t"<<sum<<"\t"<<map.size()<<"\t"<<standDev<<"\n";
        forall(const vertex &v, T){
            if(v->label != -1){
                //std::cout<<" before "<<map[v->label]<<" ";
                if(double_equals(sum,map[v->label],0.001) and min == max){
                  map[v->label] = 1;
                }
                else{                
                  map[v->label] = (map[v->label] - sum)/standDev;
                }
                //std::cout<<" after "<<map[v->label]<<"\n";
            }
        }
        //std::cout<<"\nMap size\t"<<map.size();
    }
    bool SizeCluster_SVM_Prepare::ScaleData(vvGraph &T, CellCluster &cc){  
        //Feature Scaling
       // std::cout<<"\nSizes\t"<<cc.area.size()<<"\t"<<cc.perimeter.size()<<"\t"
        //<<cc.stomataBending.size()<<"\t"<<cc.lenBr.size()<<"\t"<<cc.bending.size()<<"\t"<<cc.commonNeighbors.size()<<"\t"
        //<<cc.neighbors.size()<<"\t"<<cc.variabilityRadius.size()<<"\t"<<cc.curvature.size()<<"\t"<<cc.avgBending.size();
        FeatureScaling(T, cc.area);  
        FeatureScaling(T, cc.perimeter);  
        FeatureScaling(T, cc.bending);
        FeatureScaling(T, cc.stomataBending);
        FeatureScaling(T, cc.lenBr);
        FeatureScaling(T, cc.commonNeighbors);
        FeatureScaling(T, cc.neighbors);
        FeatureScaling(T, cc.curvature);
        FeatureScaling(T, cc.variabilityRadius);
        FeatureScaling(T, cc.stomataArea);
        FeatureScaling(T, cc.lobyConPeri);
        FeatureScaling(T, cc.lobyConAr);
        FeatureScaling(T, cc.lobyVisStom);
        FeatureScaling(T, cc.lobyVisPav);
        FeatureScaling(T, cc.lobyCirc);
        FeatureScaling(T, cc.lobyLargEmpCir);
        FeatureScaling(T, cc.avgArea);  
        FeatureScaling(T, cc.avgPerimeter);  
        FeatureScaling(T, cc.avgStomataBending);  
        FeatureScaling(T, cc.avgBending);  
        FeatureScaling(T, cc.avgCommonNeighbors);
        FeatureScaling(T, cc.avgVariabilityRadius);  
        FeatureScaling(T, cc.avgNeighbors);  
        FeatureScaling(T, cc.avgCurvature);

        //std::cout<<"\nSizes\t"<<cc.area.size()<<"\t"<<cc.perimeter.size()<<"\t"
        //<<cc.stomataBending.size()<<"\t"<<cc.lenBr.size()<<"\t"<<cc.bending.size()<<"\t"<<cc.commonNeighbors.size()<<"\t"
        //<<cc.neighbors.size()<<"\t"<<cc.variabilityRadius.size()<<"\t"<<cc.curvature.size()<<"\t"<<cc.avgBending.size();
        
        cc.FeatureMap["Area"] = cc.area;        
        cc.FeatureMap["Perimeter"] = cc.perimeter;
        cc.FeatureMap["Bending"] = cc.bending;
        cc.FeatureMap["Stomata Bending"] = cc.stomataBending;
        cc.FeatureMap["Length Breadth Ratio"] = cc.lenBr;
        cc.FeatureMap["Common Neighbors"] = cc.commonNeighbors;
        cc.FeatureMap["Neighbors"] = cc.neighbors;
        cc.FeatureMap["Curvature"] = cc.curvature;
        cc.FeatureMap["Variability Radius"] = cc.variabilityRadius;
        cc.FeatureMap["Lobyness Convex Perimeter"] = cc.lobyConPeri; 
        cc.FeatureMap["Lobyness Convex Area"] = cc.lobyConAr; 
        cc.FeatureMap["Lobyness Visibility Stomata"] = cc.lobyVisStom; 
        cc.FeatureMap["Lobyness Visibilty Pavement"] = cc.lobyVisPav; 
        cc.FeatureMap["Lobyness Circularity"] = cc.lobyCirc; 
        cc.FeatureMap["Lobyness Large Empty Circle"] = cc.lobyLargEmpCir; 
        cc.FeatureMap["Stomata Area"] = cc.stomataArea;
        cc.FeatureMap["Average Area"] = cc.avgArea;
        cc.FeatureMap["Average Perimeter"] = cc.avgPerimeter;
        cc.FeatureMap["Average Stomata Bending"] = cc.avgStomataBending;
        cc.FeatureMap["Average Bending"] = cc.avgBending;
        cc.FeatureMap["Average Common Neighbors"] = cc.avgCommonNeighbors;
        cc.FeatureMap["Average Variability Radius"] = cc.avgVariabilityRadius;
        cc.FeatureMap["Average Neighbors"] = cc.avgNeighbors;
        cc.FeatureMap["Average Curvature"] = cc.avgCurvature; 

        std::cout<<"\nFeature map size in prepare\t"<< cc.FeatureMap.size();


    }
    bool SizeCluster_SVM_Clear::run(CellCluster &cc){  
        cc.area.clear();
        cc.perimeter.clear();  
        cc.bending.clear();
        cc.stomataBending.clear();
        cc.lenBr.clear();
        cc.commonNeighbors.clear();
        cc.neighbors.clear();
        cc.curvature.clear();
        cc.variabilityRadius.clear();
        cc.stomataArea.clear();
        cc.lobyConPeri.clear();
        cc.lobyConAr.clear();
        cc.lobyVisStom.clear();
        cc.lobyVisPav.clear();
        cc.lobyCirc.clear();
        cc.lobyLargEmpCir.clear();
        cc.avgArea.clear();  
        cc.avgPerimeter.clear();  
        cc.avgStomataBending.clear();  
        cc.avgBending.clear();  
        cc.avgCommonNeighbors.clear();
        cc.avgVariabilityRadius.clear();  
        cc.avgNeighbors.clear();  
        cc.avgCurvature.clear();
        cc.SelectedFeatures.clear();
        cc.FeatureMap.clear();
        return true;
    }
    REGISTER_PROCESS(SizeCluster_SVM_Clear);

    bool SizeCluster_SVM_Prepare::run(Mesh *mesh, CellCluster &cc) //QString filename1, QString filename2, QString filename3, QString filename4,
    {  
        // Fix Corners
        FixMeshCorners FixCorners(*this);
        FixCorners.run(mesh, false, true, 0);
        // Make cells
        vvGraph T = mesh->graph();

        MeasureArea are(*this);
        are.CalculateArea(T, cc.area);

        MeasurePerimeter per(*this);
        per.CalculatePerimeter(T, cc.perimeter);

        MeasureBending bnd(*this);
        bnd.CalculateBending(T, cc.bending);

        MeasureStomatalBending stb(*this);
        stb.CalculateStomatalBending(T, cc.stomataBending); 
        
        MeasureLenBr lbr(*this);
        lbr.CalculateLenBr(T, cc.lenBr);

        MeasureNeighbors nbr(*this);
        nbr.CalculateNeighbors(T, cc.neighbors);
        std::map<int, std::set<vertex> > Nhbr;
        nbr.FindNeighbourVertices(T, Nhbr);
        
        MeasureCurvature cur(*this);
        cur.CalculateCurvature(T, cc.curvature);
        
        MeasureVariabilityRadius var(*this);
        var.CalculateVariabilityRadius(T, cc.variabilityRadius);
        
        MeasureCommonNhbrs cnb(*this);
        cnb.CalculateCommonNhbrs(T, cc.commonNeighbors);

        MeasureStomataArea sta(*this);
        sta.CalculateStomataArea(T, cc.stomataArea);

        std::cout<<"size\tcc.stomataArea\t"<<cc.stomataArea.size();

        MeasureLobyness loby(*this);
        /*loby.CalculateLobyness(T, "Convex(Perimeter)", "1", cc.lobyConPeri);
        loby.CalculateLobyness(T, "Convex(Area)", "1", cc.lobyConAr);
        loby.CalculateLobyness(T, "Visibility(Stomata)", "1", cc.lobyVisStom);
        loby.CalculateLobyness(T, "Visibility(Pavement)", "1", cc.lobyVisPav);
        loby.CalculateLobyness(T, "Circularity", "1", cc.lobyCirc);
        loby.CalculateLobyness(T, "LargestEmptySpace", "1", cc.lobyLargEmpCir);*/

        MeasureAvgArea avg(*this);
        avg.CalcAvgNhbdMeasure(T, cc.area, cc.avgArea);
        avg.CalcAvgNhbdMeasure(T, cc.perimeter, cc.avgPerimeter);
        avg.CalcAvgNhbdMeasure(T, cc.bending, cc.avgBending);
        avg.CalcAvgNhbdMeasure(T, cc.stomataBending, cc.avgStomataBending);
        avg.CalcAvgNhbdMeasure(T, cc.lenBr, cc.avgLenBr);
        avg.CalcAvgNhbdMeasure(T, cc.neighbors, cc.avgNeighbors);
        avg.CalcAvgNhbdMeasure(T, cc.curvature, cc.avgCurvature);
        avg.CalcAvgNhbdMeasure(T, cc.variabilityRadius, cc.avgVariabilityRadius);
        avg.CalcAvgNhbdMeasure(T, cc.commonNeighbors, cc.avgCommonNeighbors);
        ScaleData(T, cc);
        return true;
    }
    REGISTER_PROCESS(SizeCluster_SVM_Prepare);

    bool SizeCluster_SVM_FeatureSelection::FeatureBox(CellCluster &cc){ 
        std::cout<<"\nFeature map size"<< cc.FeatureMap.size();
        forall(TypeFeatureMap Feature, cc.FeatureMap) {
          QTreeWidgetItem *item = new QTreeWidgetItem(QStringList() << Feature.first);
          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_FeatureSelection::initialize(QStringList &parms, QWidget *parent){ 
        Mesh *m1 = currentMesh();
        m1->showHeat();
        vvGraph &T = m1->graph(); 
        QDialog dlg(parent);
        ui.setupUi(&dlg);
        FeatureBox(cc);
        std::cout<<__LINE__<<"\n";
        connect(ui.svmTreeWidget, SIGNAL(itemClicked(QTreeWidgetItem *, int)),
                this, SLOT(on_svmTreeWidget_itemClicked(QTreeWidgetItem *, int )));
        std::cout<<__LINE__<<"\n";
        connect(ui.selectAllButton, SIGNAL(clicked()), this, SLOT(on_selectAllButton_clicked()));
        std::cout<<__LINE__<<"\n";
        
        std::cout<<__LINE__<<"\n";
        dlg.exec();
        std::cout<<__LINE__<<"\n";
        selectFeatures();
        std::cout<<"Selected feaures initialize\t"<<cc.SelectedFeatures.size()<<"\n";
      }
      void SizeCluster_SVM_FeatureSelection::on_svmTreeWidget_itemClicked(QTreeWidgetItem *item, int column){
        Mesh *m1 = currentMesh();
        m1->showHeat();
        IntFloatAttr heatMap;
        float maxValue = 0, minValue = HUGE_VAL;
        vvGraph &T = m1->graph(); 
        forall(TypeFeatureMap Feature, cc.FeatureMap){
          if(item->text(0) == Feature.first){
            heatMap = Feature.second;
            break;
          }
        }
        if(heatMap.size() > 0){
          forall(to_loop_map loop, heatMap){
            if(minValue > heatMap[loop.first])
              minValue = heatMap[loop.first];
            if(maxValue < heatMap[loop.first])
              maxValue = heatMap[loop.first];
          }
          m1->wallHeat().clear();
          m1->labelHeat().clear();
          forall(const vertex v, T){
            if(v->label != -1)
              m1->labelHeat()[v->label] = heatMap[v->label];
          }
          m1->heatMapBounds() = Point2f(minValue, maxValue);
          m1->updateTriangles();
          m1->showHeat();
          m1->updateAll();
          updateState();
          updateViewer();
        }
    }
    void SizeCluster_SVM_FeatureSelection::on_selectAllButton_clicked(){
        std::cout<<__LINE__<<"\n";
        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);
        }
        std::cout<<__LINE__<<"\n";
    }
    void SizeCluster_SVM_FeatureSelection::selectFeatures(){ 
        std::vector<QString> Features;
        std::cout<<__LINE__<<"\n";
        cc.SelectedFeatures.clear();
        for(int itemCount = 0; itemCount < ui.svmTreeWidget->topLevelItemCount(); ++itemCount){
          qDebug()<<"\t"<<ui.svmTreeWidget->topLevelItem(itemCount)->text(0)<<"\t"<<ui.svmTreeWidget->topLevelItem(itemCount)->checkState(1);
          if(Qt::Checked == ui.svmTreeWidget->topLevelItem(itemCount)->checkState(1)){
              Features.push_back(ui.svmTreeWidget->topLevelItem(itemCount)->text(0));
              std::cout<<"\nFeatures size\t"<<Features.size();
          }
          
        }
        std::cout<<__LINE__<<"\n";
        forall(TypeFeatureMap Feature, cc.FeatureMap){
          forall(QString f, Features){
            if(Feature.first == f)
              cc.SelectedFeatures.push_back(Feature.second);
          }
        }
        std::cout<<__LINE__<<"\n";
        std::cout<<"Selected feaures\t"<<cc.SelectedFeatures.size()<<"\n";
    }
    REGISTER_PROCESS(SizeCluster_SVM_FeatureSelection);
    
    bool SizeCluster_SVM::run(Mesh *mesh, CellCluster &cc, QString choice)
    {  
        vvGraph &T = mesh->graph();
        if (choice == "Giant Cells"){
          forall(const vertex &v, T){
            if(v->selected and v->label != -1)
              mesh->parents()[v->label] = 1;
          }
        }
        else if (choice == "Old Stomata"){
          forall(const vertex &v, T){
            if(v->selected and v->label != -1)
              mesh->parents()[v->label] = 2;
          }
        }
        else if (choice == "Young Stomata"){
          forall(const vertex &v, T){
            if(v->selected && v->label != -1)
              mesh->parents()[v->label] = 3;
          }
        }
        else if (choice == "Large Regular Cells"){
          forall(const vertex &v, T){
            if(v->selected && v->label != -1)
              mesh->parents()[v->label] = 4;
          }
        }
        else if (choice == "Small Regular Cells"){
          forall(const vertex &v, T){
            if(v->selected && v->label != -1)
              mesh->parents()[v->label] = 5;
          }
        }
        else if (choice == "Large Marginal Cells"){
          forall(const vertex &v, T){
            if(v->selected && v->label != -1)
              mesh->parents()[v->label] = 6;
          }
        }
        else if (choice == "Small Marginal Cells"){
          forall(const vertex &v, T){
            if(v->selected && v->label != -1)
              mesh->parents()[v->label] = 7;    
          }
        }
        else if (choice == "Petiole Cells"){
          forall(const vertex &v, T){
            if(v->selected && v->label != -1)
              mesh->parents()[v->label] = 8;
          }
        }
        else if (choice == "Pavement Cells"){
          forall(const vertex &v, T){
            if(v->selected && v->label != -1)
              mesh->parents()[v->label] = 9;
          }
        }
        mesh->updateTriangles();
      return true;
    }
    REGISTER_PROCESS(SizeCluster_SVM);

    bool SizeCluster_SVM_PrepareTrainingDataFile::prepare_training_data(Mesh *mesh, CellCluster &cc, QString filename, QString choice, QString trainTest){
        std::cout<<__LINE__<<"\n";
        vvGraph &T = mesh->graph();
        if(filename.isEmpty())
          return false;
        if(!filename.endsWith(".txt", Qt::CaseInsensitive))
          filename += ".txt";
        std::cout<<__LINE__<<"\n";
        QFile file(filename);

        if(choice == "New File"){
            if(!file.open(QIODevice::WriteOnly)) 
                std::cout<<"Cannot Create a new file";
        }
        else 
            if(!file.open(QIODevice::WriteOnly | QIODevice::Append))
                std::cout<<"Cannot append to the file";
        std::cout<<__LINE__<<"\n";
        QTextStream out(&file);

        int rows = 0;
        forall(const vertex &v, T){
          int i = 0; 
          if(!(trainTest.size() and v->label > 0))
            continue;
          std::cout<<__LINE__<<"\n";
          qDebug()<<trainTest;
          std::cout<<"compare\t"<<QString::compare(trainTest, "Train", Qt::CaseSensitive);
          if(!QString::compare(trainTest, "Train", Qt::CaseSensitive)){
            if(mesh->parents()[v->label]){
              std::cout<<__LINE__<<"\n";
              rows++;
              out << mesh->parents()[v->label] <<" ";
              forall(IntFloatAttr &map, cc.SelectedFeatures){
                  forall(to_loop_map loop, map){
                    if(loop.first == v->label and !isnan(loop.second)){
                      std::cout<<__LINE__<<"\n";
                      out << i << ":"<< loop.second <<" ";
                      i++;
                    }
                  }
              }
              out << i << ":"<< "-1" <<" ";
              out << endl;               
            }
          }
          else {
            std::cout<<__LINE__<<"\n" <<"size "<<cc.SelectedFeatures.size();
              rows++;
              out << v->label <<" ";
              std::cout<<__LINE__<<"\n";
              forall(IntFloatAttr &map, cc.SelectedFeatures){
                  std::cout<<__LINE__<<"\n";
                  forall(to_loop_map loop, map){
                    std::cout<<__LINE__<<"\n";
                    if(loop.first == v->label and !isnan(loop.second)){
                      std::cout<<__LINE__<<"\n";
                      out << i << ":"<< loop.second <<" ";
                      i++;
                    }
                  }
              }
              out << i << ":"<< "-1" <<" ";
              out << endl;               
            }
        }
        file.close();        
    }
    bool SizeCluster_SVM_PrepareTrainingDataFile::run(Mesh *mesh, CellCluster &cc, QString filename, QString choice, QString trainTest){ 
        prepare_training_data(mesh, cc, filename, choice, trainTest);
    }
    REGISTER_PROCESS(SizeCluster_SVM_PrepareTrainingDataFile);

    double SizeCluster_SVM_Train::do_cross_validation()
    {
        std::cout<<__LINE__<<"\n";
        int i;
        double cross_validation;
        int total_correct = 0;
        double total_error = 0;
        double sumv = 0, sumy = 0, sumvv = 0, sumyy = 0, sumvy = 0;
        double *target = Malloc(double, prob.l);
        int nr_fold = 3;
        std::cout<<__LINE__<<"\n";
        if(&prob and &param and nr_fold and target)
        svm_cross_validation(&prob, &param, nr_fold, target);
        std::cout<<__LINE__<<"\n";
        if(param.svm_type == EPSILON_SVR ||
           param.svm_type == NU_SVR)
        {
          std::cout<<__LINE__<<"\n";
          for(i=0;i<prob.l;i++)
          {
            double y = prob.y[i];
            double v = target[i];
            total_error += (v-y)*(v-y);
            sumv += v;
            sumy += y;
            sumvv += v*v;
            sumyy += y*y;
            sumvy += v*y;
          }
          printf("Cross Validation Mean squared error = %g\n",total_error/prob.l);
          printf("Cross Validation Squared correlation coefficient = %g\n",
          ((prob.l*sumvy-sumv*sumy)*(prob.l*sumvy-sumv*sumy))/
          ((prob.l*sumvv-sumv*sumv)*(prob.l*sumyy-sumy*sumy))
          );
        }
        else
        {
          for(i=0;i<prob.l;i++){
            std::cout<<"\n"<<target[i]<<" "<<prob.y[i];
            if(target[i] == prob.y[i])
              ++total_correct;
          }
          cross_validation = (double)total_correct/double(prob.l);
          std::cout<<"Cross Validation Accuracy = \n"<<(100.0*total_correct/prob.l);
          
        }
        free(target);
                std::cout<<__LINE__<<"\n";
        return cross_validation;
    }

    bool SizeCluster_SVM_Train::run(Mesh *mesh, CellCluster &cc, QString filename, QString choice){
        //1. Prepare test data - create a test process DONE 
        //2. read the mesh into feauters & predict 

        //3. CReate training data
        //4. TRain -- fix prob.l to be dynamic DONE

        std::vector<QStringList> data;
        QFile file(filename);
        //std::cout<<readData(file, data);
        std::cout<<__LINE__<<"\n";
        std::cout<<data.size();
        std::cout<<__LINE__<<"\n";
        //std::cout<<system("svm-train -s 0 -t 0 -c 1 /home/namrata/LargeMesh/HannahData/TrainingData.txt");
        vvGraph &T = mesh->graph();
        std::cout<<__LINE__<<"\n";
        if(choice == "Train using SVM"){    
            param.svm_type = C_SVC;
            param.kernel_type = RBF; //LINEAR;
            param.degree = 3;
            param.gamma = 1/10; //0;  // 1/num_features
            param.coef0 = 0;
            param.nu = 0.5;
            param.cache_size = 100;
            param.C = 1;
            param.eps = 1e-3;
            param.p = 0.1;
            param.shrinking = 1;
            param.probability = 1;
            param.nr_weight = 0;
            param.weight_label = NULL;
            param.weight = NULL;        
            const char *error_msg;
            char *endptr;
            std::cout<<__LINE__<<"\n";
            if(!file.open(QIODevice::ReadOnly))
            {

            }
            QTextStream in(&file);
            int rowSize = 0;
            while(!in.atEnd()){
                QString line = in.readLine();
                qDebug()<<line;
                rowSize++;
            }
            file.close();
            //std::cout<<__LINE__<<"\n";
            prob.l = rowSize;
            std::cout<<"Size\n"<<prob.l;
            prob.y = Malloc(double, prob.l);
            std::cout<<__LINE__<<"\n";
            prob.x = Malloc(struct svm_node *, prob.l);
            std::cout<<__LINE__<<"\n";
            svm_node** x = Malloc(svm_node*, prob.l);
            std::cout<<__LINE__<<"\n";
            int row = 0, i = 0;
            std::cout<<__LINE__<<"\n";
            bool ok;
            if(!file.open(QIODevice::ReadOnly))
            {

            }
            //while(in.status() == QTextStream::Ok) {
            while(!in.atEnd()){
                QString line = in.readLine();
                QStringList fields = line.split(" ");  
                if(!line.size())
                    continue;
                int j = 0, flag = 0; 
                x_space = Malloc(struct svm_node, prob.l);
                prob.x[i] = &x_space[j];
                prob.y[i] = fields[0].toInt(&ok);
                std::cout<<"size"<<fields.size()<<"\n";              
                for(int p = 1; p < fields.size() - 1; p++){
                    std::cout<<__LINE__<<"\n";
                    QStringList val = fields[p].split(":");
                    //std::cout<<"\nfield"<<fields[p].toUtf8().constData()<<"\t"<<val[0].toInt(&ok)<<"\t"<<val[1].toDouble(&ok)<<"\n";
                    std::cout<<__LINE__<<"\n";
                    x_space[j].index = val[0].toInt(&ok);
                    x_space[j].value = val[1].toDouble(&ok);
                    std::cout<<"\nx_space[j]\t"<<x_space[j].index<<"\t"<<x_space[j].value;
                    std::cout<<__LINE__<<"\t j value"<<j<<"\n";
                    ++j;   
                    flag = 1;
                }
                if(flag == 1){
                    x_space[j++].index = -1;
                    std::cout<<__LINE__<<"\t j value"<<j<<"\n";
                    std::cout<<"\nx_space[j]\t"<<x_space[--j].index<<"\n";
                    i++;
                    std::cout<<__LINE__<<"\n";
                }
            }
            file.close();
            std::cout<<__LINE__<<"\n";
            error_msg = svm_check_parameter(&prob, &param);
            std::cout<<__LINE__<<"\n";
            if(error_msg)
                std::cout<<"ERROR: \n" << *error_msg;
            // Tune the model. Change C
            double cross_validation = 0;
            double Optimized_C = HUGE_VAL, max_cv = 0, Optimized_gamma = 0;
            for(double c = 1; c < 10000; c*= 4){
                for(double gamma = 0.1; gamma < 5; gamma*= 3){
                    std::cout<<"\n\n\nC: "<<c <<"\t"<<gamma;
                    param.C = c;
                    param.gamma = gamma;
                    cross_validation = do_cross_validation();
                    if(cross_validation > max_cv){
                        max_cv = cross_validation;
                        Optimized_C = param.C;
                        Optimized_gamma = param.gamma;
                    }
                }
            }
            std::cout<<"\nFinal C\t"<<param.C;
            if(Optimized_C > 0){
                param.C = Optimized_C;
                param.gamma = Optimized_gamma;
                model = svm_train(&prob, &param);
            
                std::cout<<__LINE__<<"\n";
                std::cout<<"\nFinal cross validation \t"<<max_cv<<"\t"<<Optimized_C<<"\n" <<param.gamma;
                int nr_class = svm_get_nr_class(model);
                int *labels = (int *) malloc(nr_class*sizeof(int));

                svm_get_labels(model,labels);
                printf("labels");   
                for(int p=0; p<nr_class; p++)
                    printf(" %d",labels[p]);
                printf("\n");
                free(labels);
            }
            svm_save_model((filename + ".model").toUtf8().constData(), model);
            return 0;
      }
      if(choice == "Train using Perceptron"){
         
      }
      return true;
    }
    REGISTER_PROCESS(SizeCluster_SVM_Train);
           
    bool SizeCluster_SVM_Test::run(Mesh *mesh, QString choice, CellCluster &cc)
      { 
        vvGraph &T = mesh->graph();
        std::map<int,int> count;
        //SizeCluster_SVM_Prepare pre(*this);
        //pre.run(mesh, cc);
//Problem : I need to select parameters independently from training data

        if(choice == "File"){
          std::vector<int> lab, pLab;
          std::cout<<__LINE__<<"\n";
            //std::cout<<system("svm-predict /home/namrata/MorphoDynamX_SB/mgx/branches/MorphoDynamX/AddOns/TypeRecognition/libsvm/test_data.txt /home/namrata/MorphoDynamX_SB/mgx/branches/MorphoDynamX/AddOns/TypeRecognition/libsvm/TrainingData.txt.model /home/namrata/MorphoDynamX_SB/mgx/branches/MorphoDynamX/AddOns/TypeRecognition/libsvm/Output");
            QFile file("/home/namrata/MorphoDynamX_SB/mgx/branches/MorphoDynamX/AddOns/TypeRecognition/libsvm/test_data.txt");
            
            std::cout<<__LINE__<<"\n";
            if(!file.open(QIODevice::ReadOnly))
            {
              std::cout<<__LINE__<<"\n";

            }
            QTextStream in(&file);
            bool ok;
            while(!in.atEnd()){
              std::cout<<__LINE__<<"\n";
              QString line = in.readLine();
              QStringList fields = line.split(" ");
              qDebug()<<fields[0]<<fields[0].toInt();
              lab.push_back(fields[0].toInt());
            }
            
            QFile file2("/home/namrata/MorphoDynamX_SB/mgx/branches/MorphoDynamX/AddOns/TypeRecognition/libsvm/Output");
            if(!file2.open(QIODevice::ReadOnly))
            {
              std::cout<<__LINE__<<"\n";

            }
            QTextStream in2(&file2);
            std::cout<<__LINE__<<"\n";
            while(!in2.atEnd()){
              std::cout<<__LINE__<<"\n";
              QString line = in2.readLine();
              pLab.push_back(line.toInt());
              qDebug()<< line;
              std::cout<<line.toInt()<<"\n";
            }
            file2.close();
            for(int i = 0; i<pLab.size(); i++){
                mesh->parents()[lab[i]] = pLab[i];
                std::cout<<"\n"<< lab[i] <<"\t"<< pLab[i];
            }
            std::cout<<__LINE__<<"\n";
        }
        else{
            std::cout<<__LINE__<<"\n";
            std::cout<<"Selected feauresin test\t"<<cc.SelectedFeatures.size()<<"\n";
            typedef std::pair<int,double> to_loop_IntFloatAttr;
            int svm_type = svm_get_svm_type(model);
            int nr_class = svm_get_nr_class(model);
            std::cout<<__LINE__<<"\n";
            double *prob_estimates = NULL;
            prob_estimates = (double *) malloc(nr_class*sizeof(double));
            forall(const vertex &v, T){
              std::cout<<__LINE__<<"\n";
                if(v->label == -1)
                continue;
                std::cout<<"\nlabel"<<v->label<<"\t"<<cc.SelectedFeatures.size();
                int j = 0;
                svm_node *testnode = Malloc(svm_node, prob.l);
                forall(IntFloatAttr map, cc.SelectedFeatures){
                  std::cout<<__LINE__<<"\n";
                    forall(to_loop_map loop, map){
                      std::cout<<__LINE__<<"\n";
                        if(loop.first == v->label){
                        testnode[j].index = j;
                        testnode[j].value = loop.second;
                        ++j;
                    }
                }
                testnode[j].index = -1;   //Each row of properties should be terminated with a -1 according to the readme
                double predict_label = svm_predict_probability(model,testnode,prob_estimates);    //svm_predict(model, testnode);
                std::cout<<__LINE__<<"\n";
                int estimate  = 0;
                //std::cout<<"\n\n\n\nLabel\t"<<v->label<<"\n";
                for(int p = 0;p < nr_class; p++){
                    //if(prob_estimates[p] > 0.8)
                        estimate = 1;
                    //std::cout<<predict_label<<"\t"<<prob_estimates[p] <<"\n";
                }
                if(estimate == 1)
                    mesh->parents()[v->label] = predict_label;  
                }
                std::cout<<__LINE__<<"\n";
            }
        }
        mesh->updateTriangles();
        svm_free_and_destroy_model(&model);
        svm_destroy_param(&param);
        free(prob.y);
        free(prob.x);
        free(x_space);
        return true;

    }
    REGISTER_PROCESS(SizeCluster_SVM_Test);

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

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

   
