#include <algorithm>
#include <iostream>

#include <boost/format.hpp>

#include "complexe.h"


struct special_compare : public std::unary_function<std::pair < std::pair < int, int > , std::pair < int, int > >, bool>
{
  explicit special_compare(const std::pair < std::pair < int, int > , std::pair < int, int > > &pairOfPair) :
        pairOfPair(pairOfPair) {}

  bool operator() (const std::pair < std::pair < int, int > , std::pair < int, int > > &pairOfPair2)
  {
      return (pairOfPair.first == pairOfPair2.first and pairOfPair.second == pairOfPair2.second) or (pairOfPair.first == pairOfPair2.second and pairOfPair.second == pairOfPair2.first);
  }

  std::pair < std::pair < int, int > , std::pair < int, int > > pairOfPair;
};



struct StructureComp {
  StructureComp() {}

  bool operator()(const Structure &s1, const Structure &s2) const {
          return s1.get_rna_() < s2.get_rna_();
  }
};

Complexe::Complexe(const std::vector< unsigned int > &rnaList, const std::vector< std::string > &rnaName,
                   const std::vector< std::string > &seqList,
                   const std::vector < std::pair < std::pair < unsigned int, unsigned int >, std::pair < unsigned int, unsigned int > > > &listBP,
                   std::vector < float > obj) :
    rnaList_(rnaList), rnaName_(rnaName), seqList_(seqList), listBP_(listBP), obj_(obj)
{
    std::sort(listBP_.begin(), listBP_.end());
    make_listBP2_();
}

Complexe::Complexe()
{
}

Complexe::~Complexe()
{
}

void Complexe::loadListsBP_()
{
    if( listsBP_.empty()) {
        listsBP_.push_back(listBP_);

        size_t i, j, k, size, size2, size3;
        //identify RNA with identical sequences
        std::vector< std::string > seqs;
        std::vector< std::vector< int > > idRNAs;
        //int max = 0;
        bool SYM = false;
        for(i = 0, size = seqList_.size(); i != size; i++) {
            std::vector< std::string >::iterator it = std::find(seqs.begin(), seqs.end(), seqList_[i]);
            int index = int(std::distance(seqs.begin(), it));
            if( it == seqs.end()) {
                seqs.push_back(seqList_[i]);
                std::vector<int> v { int(i) };
                idRNAs.push_back(v);
            } else {
                idRNAs[ulong(index)].push_back( int(i));
                SYM = true;
            }
        }

        if (SYM) {

            //Perform permutations for all groups of identical RNAs
            std::vector< std::vector < std::vector< int > > > permutatedRNAs;
            for(i = 0, size = idRNAs.size(); i != size; i++) {
                permutatedRNAs.push_back(std::vector < std::vector < int > > ());
                std::sort (idRNAs[i].begin(),idRNAs[i].end());

                do {
                  permutatedRNAs[i].push_back(std::vector < int > (idRNAs[i]));
                } while ( std::next_permutation(idRNAs[i].begin(),idRNAs[i].end()) );
            }

            //Mix all groups of permutated RNAs
            std::vector< std::vector < std::vector < int > > > mixedPermutations;
            std::vector< std::vector < std::vector < int > > > niveau;

            for(i = 0, size = permutatedRNAs[0].size(); i != size; i++) {
                 niveau.push_back(std::vector< std::vector < int > > {permutatedRNAs[0][i]});
            }


            std::vector< std::vector < int > > v;
            for(i = 1, size = permutatedRNAs.size(); i != size; i++) {
                mixedPermutations.clear();
                for(j = 0, size2 = permutatedRNAs[i].size(); j != size2; j++)
                    for(k = 0, size3 = niveau.size(); k != size3; k++) {
                        v = niveau[k];
                        v.push_back(permutatedRNAs[i][j]);
                        mixedPermutations.push_back(v);
                    }
                niveau = mixedPermutations;
            }

            mixedPermutations = niveau;

            //Transform into a permutation
            std::vector < std::vector < int > > permutbis, permutations;
            std::vector < int > permut;
            int max = int(rnaList_.size());
            for(i = 0, size = mixedPermutations.size(); i != size; i++) {
                permutbis = mixedPermutations[i];
                permut.clear();
                for(int j = 0; j != max; j++) {
                    k = 0;
                    while (std::find(mixedPermutations[i][k].begin(), mixedPermutations[i][k].end(), j) == mixedPermutations[i][k].end())
                        k++;
                    permut.push_back(permutbis[k][0]);
                    permutbis[k].erase(permutbis[k].begin());
                }
                permutations.push_back(permut);
            }

            //Fill listsBP_
            std::vector < std::pair < std::pair < unsigned int, unsigned int >, std::pair < unsigned int, unsigned int > > > listBPbis;
            for(i = 0, size = permutations.size(); i != size; i++) {
                for(j = 0, size2 = listBP_.size(); j != size2; j++)
                    listBPbis.push_back(std::make_pair(
                                            std::make_pair(permutations[i][listBP_[j].first.first], listBP_[j].first.second),
                                            std::make_pair(permutations[i][listBP_[j].second.first], listBP_[j].second.second) ));
                listsBP_.push_back(listBPbis);
                listBPbis.clear();
            }
        }
    }
}

std::vector< unsigned int > Complexe::get_rnaList_() const
{
    return rnaList_;
}

std::vector< std::string > Complexe::get_rnaName_() const
{
    return rnaName_;
}

std::vector< std::string > Complexe::get_seqList_() const
{
    return seqList_;
}

std::vector < float > Complexe::get_obj_() const
{
    return obj_;
}

std::vector < std::pair < std::pair < unsigned int, unsigned int > , std::pair < unsigned int, unsigned int > > > Complexe::get_listBP_() const
{
    return listBP_;
}

std::vector < std::vector < std::pair < std::pair < unsigned int, unsigned int > , std::pair < unsigned int, unsigned int > > > > Complexe::get_listsBP_() const
{
    return listsBP_;
}

void Complexe::set_obj_(std::vector < float > obj)
{
    obj_ = obj;
}

void Complexe::set_rnaList_(const std::vector< unsigned int > &rnaList)
{
    rnaList_ = rnaList;
}

void Complexe::set_rnaName_(const std::vector< std::string > &rnaName)
{
    rnaName_ = rnaName;
}

void Complexe::set_seqList_(const std::vector< std::string > &seqList)
{
    seqList_ = seqList;
}

void Complexe::make_listBP2_() {
    int rna1, rna2, i, j, i_cplx, j_cplx;
    std::tuple < int, int, int > nt1, nt2, temp;
    for (size_t k = 0, size = listBP_.size(); k != size; k++)
    {
        rna1 = int(listBP_[k].first.first);
        rna2 = int(listBP_[k].second.first);
        i = int(listBP_[k].first.second); // position of bp first nt in its rna
        j = int(listBP_[k].second.second); // position of bp second nt in its rna

        if (rna1 == 0) {
          nt1 = std::make_tuple(rna1, i, i);
        }
        else {
          i_cplx = i;
          for (size_t k = 0, size = rna1; k != size; k++) {
            i_cplx += seqList_[k].size();
          }
          nt1 = std::make_tuple(rna1, i, i_cplx);
        }

        if (rna2 == 0) {
          nt2 = std::make_tuple(rna2, j, j);
        }
        else {
          j_cplx = j;
          for (size_t k = 0, size = rna2; k != size; k++) {
            j_cplx += seqList_[k].size();
          }
          nt2 = std::make_tuple(rna2, j, j_cplx);
        }

        if (rna1 > rna2) {
            temp = nt2;
            nt2 = nt1;
            nt1 = temp;
        }

        listBP2_.push_back(std::make_pair(nt1, nt2));
    }
    std::sort(listBP2_.begin(), listBP2_.end());
}


std::string Complexe::convToDP()
{
    // sort base pairs per helices
    std::vector < std::vector < std::pair < std::tuple < int, int, int >, std::tuple < int, int, int > > > > helices;
    int rna1, rna2, i, j, i2, j2;
    bool new_hel, add;

    for (size_t k = 0, size = listBP2_.size(); k != size; k++ ) { // for bp in listBP2_ :
        i = int(std::get<2>(listBP2_[k].first)); // position of first nt in the complex
        j = int(std::get<2>(listBP2_[k].second)); // position of second nt in the complex

        new_hel = true;
        for(size_t h = 0, size2 = helices.size(); h != size2; h++) { // for hel in helices :
            add = true;
            for (size_t b = 0, size3 = helices[h].size(); b != size3; b++) { // for bp2 in hel :
                i2 = std::get<2>(helices[h][b].first);
                j2 = std::get<2>(helices[h][b].second);

                if ( (i < i2 and i2 < j and j < j2) or (i2 < i and i < j2 and j2 < j) ) {
                  add = false; // do not add bp in hel
                  break; // go to next hel
                }
            }
            if (add) {
                helices[h].push_back(listBP2_[k]); // add bp to hel
                new_hel = false; // no need to create another hel
                break; // go to next hel
            }
        }
        if (new_hel) { // bp do not fit in any helix, create a new helix
            helices.push_back(std::vector< std::pair < std::tuple < int, int, int >, std::tuple < int, int, int > > > (1, listBP2_[k]));
        }
    }

    // initialize structures to dots
    std::vector < std::string > structs;
    for (size_t k = 0, size = seqList_.size(); k != size; k++) {
        structs.push_back(std::string (ulong(seqList_[k].size()), '.'));
    }

    // choose a character for each level
    std::vector<char> openChars;
    openChars.clear();
    openChars.push_back(40); // "("
    openChars.push_back(91); // "["
    openChars.push_back(123); // "{"
    openChars.push_back(60); // "<"
    for (size_t k = 65; k != 91; k++) { // uppercase letters
        openChars.push_back(k);
    }
    openChars.push_back(92); // "\"
    for (size_t k = 48; k != 53; k++) { // "0" to "4"
        openChars.push_back(k);
    }

    int ind = 0;
    char c = openChars[ind];

    bool error = false;

    for (size_t k = 0, size = helices.size(); k != size; k++) { // hel in helices :
        for (size_t l = 0, size2 = helices[k].size(); l != size2; l++) { // bp in hel
            rna1 = std::get<0>(helices[k][l].first); // rna of first nt
            rna2 = std::get<0>(helices[k][l].second); // rna of second nt
            i = std::get<1>(helices[k][l].first); // position of first nt in its rna
            j = std::get<1>(helices[k][l].second); // position of second nt in its rna

            structs[ulong(rna1)][ulong(i)] = c; // opening char
            structs[ulong(rna2)][ulong(j)] = get_closing_char(c); // closing char
        }
        ind++;
        if ( ind == int(openChars.size()) ) {
            error = true;
            break;
        }
        else {
            c = openChars[ind];
        }
    }

    // print the new format
    std::string res = "";
    if (error) {
        res = "No dot-parenthesis display avalaible.\n";
    }
    else {
        for (size_t k = 0, size = seqList_.size(); k != size; k++) {
            res += seqList_[k] + "\n" + structs[k] + "\n";
        }
    }

    return res;
}

char Complexe::get_closing_char(char openChar)
{
    char closeChar;
    if ( openChar >= 65 and openChar <= 90 ) { // uppercase letter
        closeChar = openChar + 32; // lowercase letter
    }
    else {
        if ( openChar == 40 ) { // if openChar is a "("
            closeChar = 41; // ")"
        }
        else {
            if (openChar == 91 or openChar == 123 or openChar == 60) { // if openChar is a "[", a "{" or a "<"
                closeChar = openChar + 2;
            }
            else {
                if ( openChar == 92) { // if openChar is a "\"
                    closeChar = 47; // "/"
                }
                else { // ( openChar >= 48 and openChar <= 52 ), numbers
                    closeChar = 105 - openChar; // "0" with "9", "1" with "8", "2" with "7", "3" with "6", "4" with "5"
                }
            }
        }
    }

    return closeChar;
}


std::string Complexe::to_string()
{
    /*std::cout << "structures" << std::endl;
    for(size_t i = 0, size = strList_.size(); i != size; i++)
        std::cout << strList_[i].get_id_() << "\n" << strList_[i].convToDP() << " score = " << strList_[i].get_obj1_() << std::endl;

    std::cout << "interactions" << std::endl;
    for(size_t i = 0, size = inteList_.size(); i != size; i++)
        std::cout << inteList_[i].get_id_() << "\n" << inteList_[i].convToDP() << " score = " << inteList_[i].get_score_() << std::endl;

    std::cout << "complexe" << std::endl;*/
    return convToDP() + "\t Energy " + boost::str(boost::format("%.2f") % obj_[1])
            + "\t User constraints " + boost::str(boost::format("%.f") % obj_[0])
            //+ "\t Compatibility " + boost::str(boost::format("%.2f") % obj_[1])
            + "\t Probing data " + boost::str(boost::format("%.2f") % obj_[2]) ;
}

std::string Complexe::to_Json()
{
    return "{\"seq\":\"" + convToDP() + "\", \"energy\":\" " + boost::str(boost::format("%.2f") % obj_[1])
            + "\", \"userct\":\" " + boost::str(boost::format("%.2f") % obj_[0])
            //+ "\", \"compatibility\":\" " + boost::str(boost::format("%.2f") % obj_[1])
            + "\", \"probing\":\" " + boost::str(boost::format("%.2f") % obj_[2])+ "\"},";
}

std::string Complexe::to_forna()
{
    std::string res;
    res = "";
    for(size_t i = 0, size = seqList_.size(); i != size; i++) {
        res += ">" + rnaName_[i] + "\n" + seqList_[i] + "\n\n";
    }
    res += "[";
    for(size_t i = 0, size = listBP_.size(); i != size; i++) {
        res += "[[" + std::to_string(listBP_[i].first.first) + "," + std::to_string(listBP_[i].first.second) + "],[" + std::to_string(listBP_[i].second.first) + "," + std::to_string(listBP_[i].second.second) + "]]";
        if (i < size -1)
            res += ",";
    }
    res += "] \nEnergy " + boost::str(boost::format("%.2f") % obj_[1])
            + "\t User constraints " + boost::str(boost::format("%.f") % obj_[0])
            //+ "\t Compatibility " + boost::str(boost::format("%.2f") % obj_[1])
            + "\t Probing data " + boost::str(boost::format("%.2f") % obj_[2])
            + ";";
    return res;
}

std::string Complexe::to_graph(uint index)
{
    std::string res, res2;
    res = "t " + std::to_string(index) + "\n";
    std::vector < uint > vertices;
    uint counter = 0;
    for(size_t i = 0, size = seqList_.size(); i != size; i++) {
        for(size_t j = 0, size2 = seqList_[i].size(); j != size2; j++) {
            // Vertices
            res += "v " + std::to_string(counter) + " " + seqList_[i][j] + " " +  "\n";
            // Backbone edges
            if (j != uint(size2) - 1)
                res2 += "e " + std::to_string(counter) + " " + std::to_string(counter + 1) + " > 1\n";
            counter++;
        }
    }
    res += res2;
    // Base pair edges
    uint l1 = 0, l2 = 0;
    for(size_t i = 0, size = listBP_.size(); i != size; i++) {
        l1 = 0;
        l2 = 0;
        for(size_t j = 0, size2 = size_t(listBP_[i].first.first); j != size2; j++)
            l1 += uint(seqList_[j].size());
        for(size_t j = 0, size2 = size_t(listBP_[i].second.first); j != size2; j++)
            l2 += uint(seqList_[j].size());
        res += "e " + std::to_string(listBP_[i].first.second + l1) + " " + std::to_string(listBP_[i].second.second + l2) + " s\n";
    }
    return res;
}

bool Complexe::operator<(const Complexe &s)
{
    return (obj_[0] > s.obj_[0] //and obj_[1] > s.obj_[1]
                    and obj_[1] < s.obj_[1] and obj_[2] > s.obj_[2]);
}

bool Complexe::operator<(const Complexe &s) const
{
    return (obj_[0] > s.obj_[0] //and obj_[1] > s.obj_[1]
                    and obj_[1] < s.obj_[1] and obj_[2] > s.obj_[2]);
}

bool Complexe::operator>(const Complexe &s)
{
    return (obj_[0] < s.obj_[0] //and obj_[1] < s.obj_[1]
                    and obj_[1] > s.obj_[1] and obj_[2] < s.obj_[2]);
}

bool Complexe::operator>(const Complexe &s) const
{
    return (obj_[0] < s.obj_[0] //and obj_[1] < s.obj_[1]
                    and obj_[1] > s.obj_[1] and obj_[2] < s.obj_[2]);
}

/*bool Complexe::operator<=(const Complexe &s)
{
    return obj1_ <= s.obj1_;
}

bool Complexe::operator<=(const Complexe &s) const
{
    return obj1_ <= s.obj1_;
}

bool Complexe::operator>=(const Complexe &s)
{
    return obj1_ >= s.obj1_;
}

bool Complexe::operator>=(const Complexe &s) const
{
    return obj1_ >= s.obj1_;
}*/

bool Complexe::operator==(const Complexe& c)
{
    bool res = true;
    if( ! ( int(listBP_.size()) == int(c.listBP_.size()) and std::abs(obj_[0]-c.obj_[0]) < 0.000001
            //and std::abs(obj_[1]-c.obj_[1]) < 0.000001
            and std::abs(obj_[1]-c.obj_[1]) < 0.000001
            and std::abs(obj_[2]-c.obj_[2]) < 0.000001) ) {
        res = false;
    } else {
        for(size_t i = 0, size = listBP_.size(); i != size and res; i++)
            if (std::find(c.listBP_.begin(), c.listBP_.end(), listBP_[i]) ==  c.listBP_.end())
                res = false;
    }
    return res;
}

bool Complexe::operator==(const Complexe& c) const
{
    bool res = true;
    if( ! ( int(listBP_.size()) == int(c.listBP_.size()) and std::abs(obj_[0]-c.obj_[0]) < 0.000001
            //and std::abs(obj_[1]-c.obj_[1]) < 0.000001
            and std::abs(obj_[1]-c.obj_[1]) < 0.000001
            and std::abs(obj_[2]-c.obj_[2]) < 0.000001 ) ) {
        res = false;
    } else {
        for(size_t i = 0, size = listBP_.size(); i != size and res; i++)
            if (std::find(c.listBP_.begin(), c.listBP_.end(), listBP_[i]) ==  c.listBP_.end())
                res = false;
    }
    return res;
}

bool Complexe::operator!=(const Complexe& c)
{
    return ! ( *(this) == c );
}

bool Complexe::operator!=(const Complexe& c) const
{
    return ! ( *(this) == c );
}

std::pair < std::pair < int, int > , std::pair < int, int > > Complexe::invert(std::pair < std::pair < int, int > , std::pair < int, int > > pair) const
{
    return std::make_pair(pair.second, pair.first);
}

/*bool Complexe::isEquivalent(const Complexe &c1, const Complexe &c2)
{
    bool RES = false;
    std::vector < std::vector < std::pair < std::pair < unsigned int, unsigned int >, std::pair < unsigned int, unsigned int > > > > listsBP1 = c1.get_listsBP_();
    std::vector < std::pair < std::pair < unsigned int, unsigned int >, std::pair < unsigned int, unsigned int > > > listBP2 = c2.get_listBP_();
    size_t i = 0, j = 0, size = listsBP1.size(), size2 = listBP2.size();

   while(!RES and i != size) {
        if(listsBP1[i].size() == size2) {
            j = 0;
            do {
                RES = std::find_if(listsBP1[i].begin(), listsBP1[i].end(), special_compare(listBP2[j])) != listsBP1[i].end();
                j++;
            }
            while(RES and j != size2);
        }
        i++;
    }

    return RES;
}*/

bool Complexe::isEquivalent(const Complexe &c1, const Complexe &c2)
{
    bool RES = false;
    std::vector < std::vector < std::pair < std::pair < unsigned int, unsigned int >, std::pair < unsigned int, unsigned int > > > > listsBP1 = c1.get_listsBP_();
    std::vector < std::pair < std::pair < unsigned int, unsigned int >, std::pair < unsigned int, unsigned int > > > listBP2 = c2.get_listBP_();
    std::vector < std::pair < std::pair < unsigned int, unsigned int >, std::pair < unsigned int, unsigned int > > > listBP3;
    unsigned int rna1, rna2, locus1, locus2;
    std::pair < unsigned int, unsigned int > nt1, nt2;
    size_t i = 0, j = 0, k = 0;
    bool found;

    for (i=0; i != listsBP1.size() and RES == false; i++) {
        if (listsBP1[i].size() == listBP2.size()) {
            listBP3 = listsBP1[i];

            found = true;
            for (j=0; j != listBP2.size() and found == true; j++) { // for each bp in listBP2
                rna1 = listBP2[j].first.first;
                rna2 = listBP2[j].second.first;
                locus1 = listBP2[j].first.second;
                locus2 = listBP2[j].second.second;

                found = false;
                // check if the bp is in listsBP3
                for (k = 0; k != listBP3.size() and found == false; k++) {
                    if ( (rna1 == listBP3[k].first.first and rna2 == listBP3[k].second.first and locus1 == listBP3[k].first.second and locus2 == listBP3[k].second.second)
                          or (rna2 == listBP3[k].first.first and rna1 == listBP3[k].second.first and locus2 == listBP3[k].first.second and locus1 == listBP3[k].second.second) ) {
                        listBP3.erase(listBP3.begin()+k); // remove the bp from listBP3 so next itteration is faster
                        found  = true;
                    }
                }
            }
            if (listBP3.size() == 0) { // so if bp lists are the same
                RES = true;
            }
        }
    }

    return RES;
}