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