alg/glseqbrowser.cpp

   1 #include "alg/glseqbrowser.hpp"
   2 #include "mussa_exceptions.hpp"
   3
   4 #include <iostream>
   5 #include <stdexcept>
   6
   7 using namespace std;
   8
   9 GlSeqBrowser::GlSeqBrowser()
  10   : border_width(25),
  11     cur_ortho(400.0, 0.0, 600.0, 0.0),
  12     viewport_size(600, 400),
  13     viewport_center((cur_ortho.right-cur_ortho.left)/2+cur_ortho.left),
  14     zoom_level(2),
  15     color_mapper(),
  16     track_container()
  17 {
  18 }
  19
  20 GlSeqBrowser::GlSeqBrowser(const GlSeqBrowser& gt)
  21   : border_width(gt.border_width),
  22     max_ortho(gt.max_ortho),
  23     cur_ortho(gt.cur_ortho),
  24     viewport_size(gt.viewport_size),
  25     viewport_center(gt.viewport_center),
  26     zoom_level(gt.zoom_level),
  27     color_mapper(gt.color_mapper),
  28     track_container(gt.track_container)
  29 {
  30 }
  31
  32 void GlSeqBrowser::initializeGL()
  33 {
  34   glEnable(GL_DEPTH_TEST);
  35   glClearColor(1.0, 1.0, 1.0, 0.0);
  36   glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);
  37   glShadeModel(GL_FLAT);
  38 }
  39
  40 void GlSeqBrowser::resizeGL(int width, int height)
  41 {
  42   viewport_size.x = width;
  43   viewport_size.y = height;
  44   glViewport(0, 0, (GLsizei)width, (GLsizei)height);
  45   //update_layout();
  46 }
  47
  48 void GlSeqBrowser::paintGL() const
  49 {
  50   glClear(GL_COLOR_BUFFER_BIT|GL_DEPTH_BUFFER_BIT);
  51
  52   glPushMatrix();
  53   glMatrixMode(GL_PROJECTION);
  54   glLoadIdentity();
  55   glOrtho(cur_ortho.left, cur_ortho.right,
  56           cur_ortho.bottom, cur_ortho.top,
  57           -50.0, 50);
  58
  59   draw();
  60
  61   glPopMatrix();
  62   glFlush();
  63 }
  64
  65 void GlSeqBrowser::processSelection(GLuint hits, GLuint buffer[], GLuint bufsize)
  66 {
  67   GLuint *ptr;
  68   GLuint names;
  69   GLuint consumed_names = 0;
  70   float z1;
  71   float z2;
  72   GLuint objtype;
  73   GLuint objid;
  74   GLuint path_index = 0;
  75   GLuint pair_key_0 = 0;
  76   GLuint pair_key_1 = 0;
  77
  78   selected_paths.clear();
  79   selected_tracks.clear();
  80
  81   std::cout << "hits = " << hits << std::endl;
  82   ptr = (GLuint *) buffer;
  83   if (hits > 0)
  84     selectedMode = true;
  85   for (GLuint i=0; i < hits; ++i)
  86   {
  87     if ((i + 5) > bufsize) {
  88       std::clog << "*** selection overflow***" << std::endl;
  89     } else {
  90       consumed_names = 0;
  91       names = *ptr++;
  92       z1 = ((float)*ptr++)/0x7fffffff;
  93       z2 = ((float)*ptr++)/0x7fffffff;
  94       objtype = *ptr++; ++consumed_names;
  95       switch (objtype) {
  96         case MussaSegment:
  97           path_index = *ptr++; ++consumed_names;
  98           pair_key_0 = *ptr++; ++consumed_names;
  99           pair_key_1 = *ptr++; ++consumed_names;
 100           if (path_index < path_segments.size()) {
 101             segment_key k(pair_key_0, pair_key_1);
 102             pair_segment_map::iterator psm_i;
 103             psm_i = path_segments[path_index].find(k);
 104             if (psm_i != path_segments[path_index].end()) {
 105               Segment &seg = psm_i->second;
 106               selected_paths.insert(seg.path_ids.begin(), seg.path_ids.end());
 107             }
 108             // else something else is wrong
 109           } else {
 110             // something wasn't right
 111             clog << "invalid path_index " << path_index
 112                  << " should have been [0,"<<path_segments.size()
 113                  << ") " << endl;
 114           }
 115           break;
 116         case MussaTrack:
 117           objid = *ptr++; ++consumed_names;
 118           selected_tracks.insert(objid);
 119         break;
 120         default:
 121           cout << "unknown type " << objtype << " ";
 122           for(; consumed_names < names; ++consumed_names) {
 123             cout << consumed_names << "," << *ptr++ << " ";
 124           }
 125           cout << endl;
 126           break;
 127       }
 128     }
 129   }
 130 }
 131
 132 void GlSeqBrowser::selectRegion(int top, int left, int bottom, int right)
 133 {
 134   GLfloat x_scale = cur_ortho.width()/((float)viewport_size.x);
 135   GLfloat y_scale = cur_ortho.height()/((float)viewport_size.y);
 136   GLfloat x_left = cur_ortho.left + (left*x_scale);
 137   GLfloat x_right = right * x_scale;
 138
 139   if (top > bottom) {
 140     // woah, someone gave us a rectangle with the origin in the lower left
 141     int temp = top;
 142     bottom = top;
 143     top = temp;
 144   }
 145   // swap the orientation of canvas coordinates
 146   GLfloat y_top = cur_ortho.top-(bottom*y_scale);
 147   GLfloat y_bottom = cur_ortho.top - top * y_scale;
 148   selectedRegion = rect<float>(y_top, x_left, y_bottom, x_right);
 149
 150   // hopefully this will make a buffer big enough to receive
 151   // everything being selected
 152   //const size_t pathz_count = mussaAnalysis->paths().refined_pathz.size();
 153   //const GLuint select_buf_size = 1 + 5 * (pathz_count + sequences.size());
 154   const GLuint select_buf_size = 500000;
 155   GLuint selectBuf[select_buf_size];
 156   glSelectBuffer(select_buf_size, selectBuf);
 157   GLint hits;
 158
 159   (void)glRenderMode(GL_SELECT);
 160   glPushMatrix();
 161   glMatrixMode(GL_PROJECTION);
 162   glLoadIdentity();
 163   glOrtho(x_left, x_right, y_top, y_bottom, -50.0, 50.0);
 164   glMatrixMode(GL_MODELVIEW);
 165   glLoadIdentity();
 166
 167   draw();
 168
 169   glFlush();
 170   glPopMatrix();
 171   hits = glRenderMode(GL_RENDER);
 172   processSelection(hits, selectBuf, select_buf_size);
 173 }
 174
 175 float GlSeqBrowser::border() const
 176 {
 177   return border_width;
 178 }
 179
 180 float GlSeqBrowser::left() const
 181 {
 182   float left;
 183   if (track_container.size() == 0)
 184   {
 185     return cur_ortho.left;
 186   } else {
 187     vector<GlSequence>::const_iterator track_i = track_container.begin();
 188     left = track_i->x();
 189     for( ; track_i != track_container.end(); ++track_i)
 190     {
 191       if (track_i->x() < left)
 192         left = track_i->x();
 193     }
 194     return left-border_width;
 195   }
 196 }
 197
 198 float GlSeqBrowser::right() const
 199 {
 200   float right;
 201   if (track_container.size() == 0)
 202   {
 203     return cur_ortho.right;
 204   } else {
 205     vector<GlSequence>::const_iterator track_i = track_container.begin();
 206     right = track_i->right();
 207     for( ; track_i != track_container.end(); ++track_i)
 208     {
 209       if (track_i->right() > right)
 210         right = track_i->right();
 211     }
 212     return right+border_width;
 213   }
 214 }
 215
 216 void GlSeqBrowser::setViewportCenter(float x)
 217 {
 218   update_viewport(x, zoom_level);
 219   viewport_center = x;
 220 }
 221
 222 float GlSeqBrowser::viewportLeft() const
 223 {
 224   return cur_ortho.left;
 225 }
 226
 227 float GlSeqBrowser::viewportCenter() const
 228 {
 229   return viewport_center;
 230 }
 231
 232 float GlSeqBrowser::viewportRight() const
 233 {
 234   return cur_ortho.right;
 235 }
 236
 237 float GlSeqBrowser::viewportHeight() const
 238 {
 239   return cur_ortho.top - cur_ortho.bottom;
 240 }
 241
 242 float GlSeqBrowser::viewportWidth() const
 243 {
 244   return cur_ortho.right - cur_ortho.left;
 245 }
 246
 247 void GlSeqBrowser::setZoom(int new_zoom)
 248 {
 249   update_viewport(viewport_center, new_zoom);
 250   zoom_level = new_zoom;
 251 }
 252
 253 int GlSeqBrowser::zoom() const
 254 {
 255   return zoom_level;
 256 }
 257
 258 void GlSeqBrowser::setColorMapper(AnnotationColors& cm)
 259 {
 260   color_mapper = cm;
 261 }
 262
 263 AnnotationColors& GlSeqBrowser::colorMapper()
 264 {
 265   return color_mapper;
 266 }
 267
 268 void GlSeqBrowser::clear()
 269 {
 270   clear_links();
 271   path_segments.clear();
 272   track_container.clear();
 273 }
 274
 275 void GlSeqBrowser::push_sequence(const Sequence &s)
 276 {
 277   GlSequence gs(s, color_mapper);
 278   push_sequence(gs);
 279 }
 280
 281 void GlSeqBrowser::push_sequence(GlSequence &gs)
 282 {
 283   clear_links();
 284   pathid = 0;
 285   track_container.push_back(gs);
 286   update_layout();
 287   if (track_container.size() > 1)
 288     path_segments.push_back(pair_segment_map());
 289 }
 290
 291 const std::vector<GlSequence>& GlSeqBrowser::sequences() const
 292 {
 293   return track_container;
 294 }
 295
 296 void GlSeqBrowser::clear_links()
 297 {
 298   path_segment_map_vector::iterator psmv_i;
 299   for(psmv_i = path_segments.begin();
 300       psmv_i != path_segments.end();
 301       ++psmv_i)
 302   {
 303     psmv_i->clear();
 304   }
 305 }
 306
 307 void
 308 GlSeqBrowser::link(const vector<int>& path, const vector<bool>& rc, int length)
 309 {
 310   if (path.size() < 2) {
 311     // should i throw an error instead?
 312     return;
 313   }
 314   if (path.size() != rc.size()) {
 315     throw runtime_error("path and reverse compliment must be the same length");
 316   }
 317   vector<int>::const_iterator path_i = path.begin();
 318   vector<bool>::const_iterator rc_i = rc.begin();
 319   int track_i = 0;
 320   int prev_x = *path_i; ++path_i;
 321   bool prev_rc = *rc_i; ++rc_i;
 322   while (path_i != path.end() and rc_i != rc.end())
 323   {
 324     segment_key p(prev_x, *path_i);
 325     pair_segment_map::iterator found_segment = path_segments[track_i].find(p);
 326     if (found_segment == path_segments[track_i].end()) {
 327       // not already found
 328       float y1 = track_container[track_i].y();
 329             y1 -= track_container[track_i].height()/2;
 330       float y2 = track_container[track_i+1].y();
 331             y2 += track_container[track_i+1].height()/2;
 332
 333       Segment s(prev_x, y1, *path_i, y2, prev_rc);
 334       s.path_ids.insert(pathid);
 335       path_segments[track_i][p] = s;
 336     } else {
 337       //found
 338       found_segment->second.path_ids.insert(pathid);
 339     }
 340     prev_x = *path_i;
 341     prev_rc = *rc_i;
 342     ++track_i;
 343     ++path_i;
 344     ++rc_i;
 345   }
 346   // pathid is reset by push_sequence
 347   ++pathid;
 348 }
 349
 350 const set<int>& GlSeqBrowser::selectedPaths() const
 351 {
 352   return selected_paths;
 353 }
 354
 355 void GlSeqBrowser::centerOnPath(const vector<int>& paths)
 356 {
 357   if (paths.size() != track_container.size()) {
 358     throw mussa_error("Path length didn't match the number of sequences");
 359   }
 360
 361   for(size_t track_i = 0; track_i != track_container.size(); ++track_i)
 362   {
 363     track_container[track_i].setX(viewport_center - paths[track_i]);
 364   }
 365 }
 366
 367 void GlSeqBrowser::update_viewport(float center, int new_zoom)
 368 {
 369   float max_width = max_ortho.width();
 370   // division by zero is a major bummer
 371   if (new_zoom < 1) {
 372     new_zoom = 1;
 373   }
 374   float new_max_width = max_width / new_zoom;
 375   cur_ortho.left = center-new_max_width;
 376   cur_ortho.right = center+new_max_width;
 377 }
 378
 379 void GlSeqBrowser::update_layout()
 380 {
 381   typedef std::vector<GlSequence>::iterator glseq_itor_type;
 382   float available_height = (float)cur_ortho.top - 2 * (float)border_width;
 383   float max_base_pairs = 0;
 384   size_t track_count = track_container.size();
 385
 386   if (track_count > 1) {
 387     // we have several sequences
 388     float track_spacing = available_height / (track_count-1);
 389     float y = available_height + (float)border_width;
 390     for(glseq_itor_type seq_i = track_container.begin();
 391         seq_i != track_container.end();
 392         ++seq_i, y-=track_spacing)
 393     {
 394       seq_i->setX(0);
 395       seq_i->setY(y);
 396       if (seq_i->length() > max_base_pairs)
 397         max_base_pairs = seq_i->length();
 398     }
 399   } else if (track_count == 1) {
 400     // center the single track
 401     glseq_itor_type seq_i = track_container.begin();
 402     seq_i->setX(0);
 403     seq_i->setY(viewport_size.x /2);
 404     max_base_pairs = seq_i->length();
 405   } else {
 406     // nothing to do as we're empty
 407     return;
 408   }
 409   max_ortho.right = max_base_pairs + border_width;
 410   max_ortho.left = -border_width;
 411   max_ortho.top = viewport_size.x;
 412   max_ortho.bottom = 0;
 413   cur_ortho = max_ortho;
 414   viewport_center = (cur_ortho.width()/2) + cur_ortho.left;
 415 }
 416
 417 void GlSeqBrowser::draw() const
 418 {
 419   glMatrixMode(GL_MODELVIEW);
 420   glInitNames();
 421   glPushName(MussaSegment);
 422   draw_segments();
 423   glLoadName(MussaTrack);
 424   draw_tracks();
 425   glPopName();
 426   // a selection shouldn't have a glName associated with it
 427   draw_selection();
 428 }
 429
 430 void GlSeqBrowser::draw_selection() const
 431 {
 432   // draw selection box
 433   glEnable(GL_BLEND);
 434   glDepthMask(GL_FALSE);
 435   if (selectedMode) {
 436     glColor4f(0.6, 0.6, 0.6, 0.9);
 437     glRectf(selectedRegion.left, selectedRegion.top,
 438             selectedRegion.right, selectedRegion.bottom);
 439   }
 440   glDepthMask(GL_TRUE);
 441   glDisable(GL_BLEND);
 442 }
 443
 444 void GlSeqBrowser::draw_tracks() const
 445 {
 446   for(size_t track_i = 0; track_i != track_container.size(); ++track_i)
 447   {
 448     glPushName(track_i);
 449     track_container[track_i].draw(cur_ortho.left, cur_ortho.right);
 450     glPopName();
 451   }
 452 }
 453
 454 void GlSeqBrowser::draw_segments() const
 455 {
 456   glLineWidth(0.1);
 457   // each vector contains path_segment_maps of all the connections
 458   // between this track and the next
 459   path_segment_map_vector::const_iterator psmv_i;
 460   for(psmv_i = path_segments.begin();
 461       psmv_i != path_segments.end();
 462       ++psmv_i)
 463   {
 464     path_segment_map_vector::difference_type path_index;
 465     path_index = psmv_i - path_segments.begin();
 466     // these maps contain the pair index (used so we dont keep drawing the
 467     // same segment) and the actual segment structure.
 468     pair_segment_map::const_iterator psm_i;
 469     for(psm_i = psmv_i->begin();
 470         psm_i != psmv_i->end();
 471         ++psm_i)
 472     {
 473       // grab the index into our segment map
 474       const segment_key& key = psm_i->first;
 475       // the second element of our map pair is a segment
 476       const Segment &s = psm_i->second;
 477       // need to do something so we can detect our selection
 478       vector<int> selected;
 479       set_intersection(selected_paths.begin(), selected_paths.end(),
 480                        s.path_ids.begin(), s.path_ids.end(),
 481                        back_inserter(selected));
 482
 483       if (not s.reversed) {
 484         if (selected_paths.size() == 0 or selected.size() > 0) {
 485           glColor3f(1.0, 0.0, 0.0);
 486         } else {
 487           glColor3f(1.0, 0.8, 0.8);
 488         }
 489       } else {
 490         if (selected_paths.size() == 0 or selected.size() > 0) {
 491           glColor3f(0.0, 0.0, 1.0);
 492         } else {
 493           glColor3f(0.8, 0.8, 1.0);
 494         }
 495       }
 496       // save the multipart name for our segment
 497       glPushName(path_index); glPushName(key.first); glPushName(key.second);
 498       glBegin(GL_LINES);
 499       glVertex3f(s.start.x, s.start.y, -1);
 500       glVertex3f(s.end.x, s.end.y, -1);
 501       glEnd();
 502       // clear the names
 503       glPopName(); glPopName(); glPopName();
 504     }
 505   }
 506 }