//Kevin Lynagh
//October 2009
//Use the HTML5 Canvas element to draw matrix plots

//global variables
var n; //how big is our square matrix?
var cms = new Array; //an array of binary matricies with indicies that correspond to the table rows. When the cutoff sliders change, update the matricies in here -- these are the ones that are plotted.
var pes = new Array; //likewise for the principal eigenvektors
var epsilon = 0.01; // if the pe changes less than this on an iteration of the power method, consider it close enough
var dots = new Array;
var px = 4; //how big should the matrix plot be?

//keep all of the contact maps in here, with the following fixed names corresponding to the offered flavors. It feels a bit gross to hardcode these, but trying to abstract them might be a bit too tricky
var cm_hashes = new Array;
var pe_hashes = new Array; //keep track of the original PEs returned by the server in case we accidently zero the one in the pes array.
//radial_alpha_carbons
//radial_sidechain_centres
//atomwise_voronoi
//radial measures have all the pairwise distances A, voronoi have interfacial area (A^2)

//vars for browser object caching
var pe_display, matrix_display, pdb_id_select, dot_display, cm_table;

//what pdb are we looking for?
var pdb_id;
var seqnos;
var missing_residues;


function cm_init(){

  //add a change event handler for the initial row
  $($('#cm_table .cm_row:last').find('.cm_type_select')).change(function(){update_cm_table_row(0);});
  $($('#cm_table .cm_row:last').find('.mathematica_link')).click(function(){row_mathematica_popup(0);});

  //cache objects to avoid slow DOM lookups in the browser
  dot_display = $('#dot_display');
  cm_table = $('#cm_table');
  pdb_id_select = $('#pdb_id_select');
  pdb_id = pdb_id_select.val();
  pe_display = $('#pe')[0];
  matrix_display = $('#matrix')[0];


  $("#pixel_slider").slider({animate: true,
                             min: 1,
                             max: 10,
                             value:4,
                             slide: function(e, ui){
                               px = ui.value;
                               redraw_plots();
                             }
                            });
  show_pdb(pdb_id_select.val());
}

//change the displayed pdb
function show_pdb(pdb){
  //reset global vars
  cms = new Array;
  pes = new Array;
  cm_hashes = new Array;
  pe_hashes = new Array;

  pdb_id = pdb;
  $('#cm_table').find("tr:gt(1)").remove();


  //get general info about this pdb, start loading up table elements
  $.getJSON( '/json/' + pdb_id, function(p){
              seqnos = p.seqnos;
              missing_residues = p.missing_residues;
              n = seqnos.length;
              update_cm_table_row(0);
            });
}


//each row in the cm_table corresponds to matrix displays and PE display canvas elements. When we add or delete a row, add or remove those elements from the page. I suspect it will be much faster (in terms of slider-response feel) to have each cm on its own canvas. Though for simplicity right now, I'll just have one element
function new_cm_table_row(){
  var tr = cm_table.find('.cm_row:last').clone();
  $('#cm_table .cm_row:last').after(tr);
  var new_row_index = $('#cm_table .cm_row').length-1;
  //add a change handler
  $($('#cm_table .cm_row:last').find('.cm_type_select')).change(function(){update_cm_table_row(new_row_index);});
  $($('#cm_table .cm_row:last').find('.mathematica_link')).click(function(){row_mathematica_popup(new_row_index);});

  //remind this new (well all of, actually) color picker
  new jscolor.bind;
  update_cm_table_row(new_row_index);
}

function update_cm_table_row(row_index){
  var row = $('#cm_table .cm_row')[row_index];
  var cutoff_input = $($(row).find('.cutoff')[0]);
  var cm_type = $($(row).find('.cm_type_select option:selected')[0]).val();


  //destroy the old slider,
  $($(row).find('.cutoff_slider')[0]).slider('destroy');
  $($(row).find('.cutoff_slider')[0]).html("<div class='cutobff_slider'></div>");

  // Create slider defaults based on the row type and set the cutoff_input
  var new_slider = $($(row).find('.cutoff_slider')[0]);
  switch(cm_type){
  case 'radial_sidechain_centres': //same parameters as radial_alpha_carbons, let this case fall through
  case 'radial_alpha_carbons':
    new_slider.slider({animate: true,
                       range: 'min',
                       min: 40,
                       max: 250,
                       value:70,
                       slide: function(e, ui){
                         value = ui.value;
                         // update the input box and redraw plots
                         cutoff_input.val( value/10.0 );
                         redraw_plots();
                       }
                      });
    cutoff_input.val(7.0);
    break;

  case 'atomwise_voronoi':
    new_slider.slider({animate: true,
                       range: 'min',
                       min: 0,
                       max: 200,
                       value:10,
                       slide: function(e, ui){
                         value = ui.value;
                         // update the input box and redraw plots
                         cutoff_input.val(value);
                         redraw_plots();
                       }});
    cutoff_input.val(10);
    break;
  }
  $(new_slider).slider('disable');

  //get data if we don't have it cached already
  if (cm_hashes[cm_type] == undefined){
    $.getJSON( '/json/' + pdb_id + '/' + cm_type, function(j){
                cm_hashes[cm_type] = j.cm;
                pe_hashes[cm_type] = j.pe;
                pes[row_index] = jQuery.extend(true, [], j.pe);
                redraw_plots();
                //the slider starts off as disabled until we get the data
                $(new_slider).slider('enable');
              });
  }else{
    redraw_plots();
    $(new_slider).slider('enable');
  }
}





//this gets called when the cutoff distance or pixel size changes, and redraws all the canvases
function redraw_plots(){

  //make the canvas elements larger if we need to..
  pe_display.width = n*px;
  matrix_display.width =  n*px;
  matrix_display.height = n*px;

  //draw
  clear_canvases();
  $(cm_table).find('.cm_row').each(function(i, row){
                                     var cutoff = $(row).find('.cutoff')[0].value;
                                     var cm_type = $(row).find('.cm_type_select')[0].value;
                                     var rgb_color = $(row).find('.color')[0].color.rgb;
                                     //the colors from the picker are between 0 and 1, scale them to 0--255
                                     rgb_color = $.map(rgb_color, function(x){return Math.floor(x*255);});
                                     var color = 'rgba(' + rgb_color.join(',') + ', 0.6)';

                                     //if the cm_type is voronoi, we want the cutoff to exclude values LESS than the facial area, otherwise we want the cutoff to get rid of values with MORE than the value (the distance)
                                     if (cm_type == 'atomwise_voronoi'){
                                       var cutoff_function = function(x){
                                         if(x > cutoff){
                                           return 1;
                                         }else{
                                           return 0;
                                         }
                                       };
                                     }else{
                                       var cutoff_function = function(x){
                                         if(x < cutoff){
                                           return 1;
                                         }else{
                                           return 0;
                                         }
                                       };
                                     }

                                     cms[i] = get_cm(cm_type, cutoff_function);
                                     pes[i] = get_pe(cms[i]);
                                     //console.log(pes[i].slice(0, 10));
                                     draw_matrix(cms[i], color);
                                     draw_pe(pes[i], color);

                                     //these are shown on table hover
                                     //dots[i] = dot(pe, vpe);
                                   });

}

function get_cm(cm_type, cutoff_function){
  //  console.log('calculating ' + cm_type + ' ' + cutoff);
  // Deep copy from the hash to create a new cm
  var cm = jQuery.extend(true, [], cm_hashes[cm_type]);
  for (i=0; i<n; i++){
    for (j=0; j<n; j++){
      //if something is missing, it will have a -1
      if (cm[i][j] == -1){
        cm[i][j] = 0;
      }else{
        cm[i][j] = cutoff_function(cm[i][j]);
      }
    }
  }
  return cm;
}




function clear_canvases(){
  var elem = matrix_display;
  if (elem && elem.getContext) {
    var context = elem.getContext('2d');
    if (context)
      context.clearRect(0, 0, elem.width, elem.height);
  }
  elem = pe_display;
  if (elem && elem.getContext) {
    var context = elem.getContext('2d');
    if (context)
      context.clearRect(0, 0, elem.width, elem.height);
    //flip the coordinate system so positive y values go upwards
    context.scale(1, -1);
    context.translate(0, -elem.height);

    //scale to make PE components taller
    context.scale(1, 200);
  }
}

function draw_matrix(m, color){
  var elem = matrix_display;
  if (elem && elem.getContext) {
    var context = elem.getContext('2d');
    if (context) { // okay, lets start drawing.
      context.fillStyle = color;
      //loop through the contact map and fill a pixel with the color when there is a one at that position in the matrix
      $(m).each(function(i, row){
                  //go through all of the values on the row
                  $(row).each(function(j, val){
                                if (val == 1){
                                  //x,y (y is defined to point downwards, width, height)
                                  context.fillRect(i*px, j*px, px, px);
                                }
                              });
                }
               );
    }
  }
}


function draw_pe(pe, color){
  var elem = pe_display;
  if (elem && elem.getContext) {
    var context = elem.getContext('2d');
    if (context) {
      //note that this drawing context is flipped upsidedown and scaled in clear_canvases()
      context.fillStyle = color;

      //draw a bar for each radial PE
      $(pe).each(function(i, pe_component){
                   context.fillRect(i*px, 0, px, pe_component);
                 });
    }
  }
}










//output a contact matrix as a mathematica file
function matrix_to_mathematica(matrix, comment){
  var str = "(* Kevin Lynagh; dirigibleFlightcraft.com *)\n";
  str += "(* use N[cm] to have Mathematica use much faster floating point maths. Also try ArrayPlot[cm] for a picture and Eigensystem[cm] for the spectrum\n";
  if(comment != undefined)
    str += '(* ' + comment + " *)\n\n";
  str += 'cm = {' +
    jQuery.map(matrix, function(row, i){
                 return '{' + row.join(', ') + '}';
               }).join(",\n") +
    '};';

  return str;
}

function row_mathematica_popup(row_index){
  var row = $('#cm_table .cm_row')[row_index];
  var cutoff = $($(row).find('.cutoff')[0]).val();
  var cm_type = $($(row).find('.cm_type_select option:selected')[0]).val();

  var comment;
  switch(cm_type){
  case 'radial_sidechain_centres':
    comment = '' + cutoff + ' angstrom cutoff for radial contact map of sidechain centroids'; break;
  case 'radial_alpha_carbons':
    comment = '' + cutoff + ' angstrom cutoff for radial contact map of alpha carbons'; break;
  case 'atomwise_voronoi':
    comment = '' + cutoff + ' angstrom squared interfacial area cutoff for atomwise Voronoi contact map'; break;
  }
  comment += ' of pdb ' + pdb_id;
  var str = matrix_to_mathematica(cms[row_index], comment);
  popup(str);
}


function popup(term)  // write corresponding content to the popup window
{
  var popupWin = window.open("", "puWin",  "width=480,height=200,scrollbars,dependent,resizable");
  popupWin.document.open("text/plain");
  popupWin.document.write(term);
  popupWin.document.close();  // close layout stream
  popupWin.focus();  // bring the popup window to the front
}































//////////////////////
//Linear Algebra
//////////////////////

function magnitude(a){
  var grr = 0;
  n = a.length;
  for(var i=0; i<n; i++)
    grr += a[i]*a[i];

  return Math.sqrt(grr);
}

//matrix multiply
function left_multiply(matrix, vector){
  var new_vector = new Array(n);
  for(var i=0; i<n; i++){
    new_component = 0;
    for(var j=0; j<n; j++){
      new_component += matrix[i][j] * vector[j];
    }
    new_vector[i] = new_component;
  }
  return new_vector;
}

function normalize(vector){
  var normalized_vector = new Array(n);
  mag = magnitude(vector);
  for(var i=0; i<n; i++){
    normalized_vector[i] = vector[i] / mag;
  }
  return normalized_vector;
}

function dot(v1, v2){
  var product = 0;
  n = v1.length;
  for(var i=0; i<n; i++){
    product += v1[i]*v2[i];
  }
  return product;
}

//use the power method to find the principal eigenvector of a matrix
function get_pe(matrix, seed_vector){
  if (seed_vector == undefined) {
    var seed_vector = new Array(n);
    for(var i=0; i<n; i++)
      //set the seed vector elements to unity if we aren't given a better guess.
      seed_vector[i] = 1;
  }
  seed_vector = normalize(seed_vector);

  //for each iteration of the algorithm, just multiply the old vector again
  for(var l=0; l<50; l++){
    new_vector = left_multiply(matrix, seed_vector);

    //make sure we don't accidently erase the radial PE
    // if (magnitude(new_vector) == 0){
    //   alert('Contact distance too small, matrix is 0. Reverting to 7 angstrom eigenvektor');
    //   $('#input_radial_cutoff').val(7);
    //   pe = orig_pe;
    //   redraw_plots();
    //   return pe;
    // }
    new_vector = normalize(new_vector);


    //break if the new vector and old vector are converging (dot product differs from unity less than some epsilon)
    if(1-dot(seed_vector, new_vector) < epsilon){
      //console.log('converged after ' + l + ' iterations');
      break;
    }
    seed_vector = new_vector;
  }
  return new_vector;
}