qui/MussaAlignedWindow.cpp

   1 #include <list>
   2 #include <vector>
   3 #include <sstream>
   4
   5 #include <QVBoxLayout>
   6
   7 #include "qui/MussaAlignedWindow.hpp"
   8 #include "alg/sequence.hpp"
   9
  10 using namespace std;
  11
  12 MussaAlignedWindow::MussaAlignedWindow(Mussa& m,
  13                                        const set<int>& sel_paths,
  14                                        QWidget *parent)
  15   : QWidget(parent),
  16     analysis(m)//,
  17     //selected_paths(sel_paths)
  18
  19 {
  20   browser.setSequences(analysis.sequences(), analysis.colorMapper());
  21   setSelectedPaths(m, sel_paths);
  22   setAlignment(0);
  23   browser.zoomToSequence();
  24   computeMatchLines();
  25
  26   QVBoxLayout *layout = new QVBoxLayout;
  27   layout->addWidget(&browser);
  28   layout->addWidget(&status);
  29   setLayout(layout);
  30
  31   ostringstream message;
  32   message << "Selected " << selected_paths.size() << " paths";
  33   status.showMessage(message.str().c_str(), 5000);
  34 }
  35
  36
  37 void MussaAlignedWindow::setSelectedPaths(Mussa &m, const set<int>& sel_paths)
  38 {
  39   // sets are sorted
  40   set<int>::iterator sel_i = sel_paths.begin();
  41   list<ExtendedConservedPath>::const_iterator path_i = m.paths().rpbegin();
  42   list<ExtendedConservedPath>::const_iterator path_end = m.paths().rpend();
  43   size_t path_size = m.paths().refined_pathz.size();
  44   size_t pathid=0;
  45
  46   selected_paths.reserve(sel_paths.size());
  47   while (pathid != path_size and sel_i != sel_paths.end())
  48   {
  49     assert (*sel_i >= 0);
  50     size_t sel_pathid = (size_t)(*sel_i);
  51     if (pathid == sel_pathid) {
  52       cout << "found one " << *path_i << endl;
  53       selected_paths.push_back(*path_i);
  54       ++pathid;
  55       ++path_i;
  56       ++sel_i;
  57     } else if (pathid < sel_pathid) {
  58       ++pathid;
  59       ++path_i;
  60     } else if (pathid > sel_pathid) {
  61       ++sel_i;
  62     }
  63   }
  64 }
  65
  66 void MussaAlignedWindow::setAlignment(size_t alignment_index)
  67 {
  68   if (selected_paths.size() > 0) {
  69     browser.centerOnPath(selected_paths[alignment_index].track_indexes);
  70   }
  71 }
  72
  73 void MussaAlignedWindow::computeMatchLines()
  74 {
  75   const vector<Sequence>& raw_sequence = analysis.sequences();
  76   vector<int> aligned_path;
  77   int i2, i3;
  78   int x_start, x_end;
  79   int window_length, win_i;
  80   int rc_1 = 0;
  81   int rc_2 = 0;
  82   bool rc_color;
  83   vector<bool> rc_list;
  84   bool full_match;
  85   vector<bool> matched;
  86   int align_counter;
  87
  88   align_counter = 0;
  89   for(vector<ExtendedConservedPath >::iterator pathz_i=selected_paths.begin();
  90       pathz_i != selected_paths.end();
  91       ++pathz_i)
  92   {
  93     if (true /* show_aligns[align_counter] */)
  94     {
  95       ExtendedConservedPath& a_path = *pathz_i;
  96       window_length = a_path.window_size;
  97       // determine which parts of the path are RC relative to first species
  98       rc_list = a_path.reverseComplimented();
  99
 100       // loop over each bp in the conserved region for all sequences
 101       for(win_i = 0; win_i < window_length; win_i++)
 102       {
 103         aligned_path.clear();
 104         // determine which exact base pairs match between the sequences
 105         full_match = true;
 106         for(i2 = 0; i2 < a_path.size()-1; i2++)
 107         {
 108           // assume not rc as most likely, adjust below
 109           rc_1 = 0;
 110           rc_2 = 0;
 111           // no matter the case, any RC node needs adjustments
 112           if (a_path[i2] < 0)
 113             rc_1 = window_length-1;
 114           if (a_path[i2+1] < 0)
 115             rc_2 = window_length-1;
 116
 117           x_start = (abs(a_path[i2]-rc_1+win_i));
 118           x_end =   (abs(a_path[i2+1]-rc_2+win_i));
 119
 120           // RC case handling
 121           // ugh, and xor...only want rc coloring if just one of the nodes is rc
 122           // if both nodes are rc, then they are 'normal' relative to each other
 123           if ( ( rc_list[i2] || rc_list[i2+1] ) &&
 124               !(rc_list[i2] && rc_list[i2+1] ) )
 125           { //the hideous rc matching logic - not complex, but annoying
 126             if ( !( ( (raw_sequence[i2][x_start] == 'A') &&
 127                     (raw_sequence[i2+1][x_end] == 'T') ) ||
 128                   ( (raw_sequence[i2][x_start] == 'T') &&
 129                     (raw_sequence[i2+1][x_end] == 'A') ) ||
 130                   ( (raw_sequence[i2][x_start] == 'G') &&
 131                     (raw_sequence[i2+1][x_end] == 'C') ) ||
 132                   ( (raw_sequence[i2][x_start] == 'C') &&
 133                     (raw_sequence[i2+1][x_end] == 'G') ) ) )
 134               full_match = false;
 135           }
 136           else
 137           {
 138             if (!( (raw_sequence[i2][x_start] == raw_sequence[i2+1][x_end]) &&
 139                   (raw_sequence[i2][x_start] != 'N') &&
 140                   (raw_sequence[i2+1][x_end] != 'N') ) )
 141               full_match = false;
 142           }
 143         }
 144
 145         // draw for matches stretching across all sequences
 146         if (full_match)
 147         {
 148           //fl_line_style(FL_SOLID, 1, NULL);
 149
 150           // now can draw the line for each bp in this window that matches
 151           // grrr, need to ask if anyone cares if I switch the seq
 152           // top-bot order...
 153           i3 = 0;
 154           //y_loc = y_min + 5;
 155           for(i2 = 0; i2 < a_path.size()-1; i2++)
 156           {
 157             // assume not rc as most likely, adjust below
 158             rc_1 = 0;
 159             rc_2 = 0;
 160             // no matter the case, any RC node needs adjustments
 161             if (a_path[i2] < 0)
 162             {
 163               rc_1 = window_length;
 164             }
 165             if (a_path[i2] < 0)
 166             {
 167               rc_2 = window_length;
 168             }
 169
 170             // set offset based on current alignment for which bp to show
 171             //offset1 = seq_align_offsets[i2];
 172             //offset2 = seq_align_offsets[i2+1];
 173
 174             if ( ( rc_list[i2] || rc_list[i2+1] ) &&
 175                 !(rc_list[i2] && rc_list[i2+1] ) ) {
 176               ;
 177               //fl_color(FL_BLUE);
 178             } else {
 179               ;
 180               //fl_color(FL_RED);
 181             }
 182
 183             // maybe shouldn't recalc these, but store values from first loop
 184             x_start = (abs((int) (a_path[i2]-rc_1+win_i)));
 185             x_end =   (abs((int) (a_path[i2+1]-rc_2+win_i)));
 186             aligned_path.push_back(x_start);
 187             // if we're on the last time through the loop, save x_end too
 188             if (i2 == a_path.size()-2) {
 189               aligned_path.push_back(x_end);
 190             }
 191           }
 192         }
 193         if (aligned_path.size() > 0) {
 194           browser.link(aligned_path, rc_list,1);
 195         }
 196       }
 197     }
 198     align_counter++;
 199   }
 200 }
 201