src/jalview/datamodel/CigarBase.java

   1 /*
   2  * Jalview - A Sequence Alignment Editor and Viewer (Version 2.8)
   3  * Copyright (C) 2012 J Procter, AM Waterhouse, LM Lui, J Engelhardt, G Barton, M Clamp, S Searle
   4  *
   5  * This file is part of Jalview.
   6  *
   7  * Jalview is free software: you can redistribute it and/or
   8  * modify it under the terms of the GNU General Public License
   9  * as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.
  10  *
  11  * Jalview is distributed in the hope that it will be useful, but
  12  * WITHOUT ANY WARRANTY; without even the implied warranty
  13  * of MERCHANTABILITY or FITNESS FOR A PARTICULAR
  14  * PURPOSE.  See the GNU General Public License for more details.
  15  *
  16  * You should have received a copy of the GNU General Public License along with Jalview.  If not, see <http://www.gnu.org/licenses/>.
  17  */
  18 package jalview.datamodel;
  19
  20 import java.util.*;
  21
  22 public abstract class CigarBase
  23 {
  24   /**
  25    * Base class for compact idiosyncratic representation of gaps and aligned
  26    * residues Regards to Tom Oldfield for his DynamicArray class. 17th July 2006
  27    * Not thread safe.
  28    */
  29   public CigarBase()
  30   {
  31     // nothing to be done (probably)
  32   }
  33
  34   protected int length = 0;
  35
  36   protected int _inc_length = 10; // extension range for addition of new
  37
  38   // operations
  39
  40   protected char[] operation = null;
  41
  42   protected int[] range = null;
  43
  44   /**
  45    * Range of Hidden residues in seq (translated as deleted in seq)
  46    */
  47   public static final char D = 'D';
  48
  49   /**
  50    * Range of insertions to seq
  51    */
  52   public static final char I = 'I';
  53
  54   /**
  55    * Range of aligned residues
  56    */
  57   public static final char M = 'M';
  58
  59   static protected final char _case_shift = 'a' - 'A';
  60
  61   /**
  62    * Ugly function to get edited sequence string, start and end symbol positions
  63    * and the deletion regions as an array of int pairs May return null for an
  64    * empty cigar string. May return null for deletion ranges if there are none.
  65    *
  66    * @param reference
  67    *          - the symbol sequence to apply the cigar operations to (or null if
  68    *          no sequence)
  69    * @param GapChar
  70    *          - the symbol to use for Insert operations
  71    * @return Object[] { String, int[] {start, startcol, end, endcol}, int[][3]
  72    *         {start, end, col} or null} the gapped sequence, first and last
  73    *         residue index, and the deletion ranges on the reference sequence
  74    */
  75   public Object[] getSequenceAndDeletions(String reference, char GapChar)
  76   {
  77     int rlength = 0;
  78     int[][] deletions = new int[length][];
  79     int[][] trunc_deletions = null;
  80     StringBuffer sq = new StringBuffer();
  81     int cursor = 0, alcursor = 0, start = 0, startpos = 0, end = 0, endpos = 0, delcount = -1;
  82     boolean consecutive_del = false;
  83     if (length == 0)
  84     {
  85       return null;
  86     }
  87     if (reference != null)
  88     {
  89       rlength = reference.length();
  90     }
  91     boolean modstart = true;
  92     for (int i = 0; i < length; i++)
  93     {
  94       switch (operation[i])
  95       {
  96       case D:
  97         if (!consecutive_del)
  98         {
  99           deletions[++delcount] = new int[]
 100           { cursor, 0, alcursor };
 101         }
 102         cursor += range[i];
 103         deletions[delcount][1] = cursor - 1;
 104         consecutive_del = true;
 105         break;
 106       case I:
 107         consecutive_del = false;
 108         for (int r = 0; r < range[i]; r++)
 109         {
 110           sq.append(GapChar);
 111           alcursor++;
 112         }
 113         break;
 114       case M:
 115         consecutive_del = false;
 116         if (modstart)
 117         {
 118           start = cursor;
 119           startpos = alcursor;
 120           modstart = false;
 121         }
 122         if (reference != null)
 123         {
 124           int sbend = cursor + range[i];
 125           if (sbend > rlength)
 126           {
 127             sq.append(reference.substring(cursor, rlength));
 128             while (sbend-- >= rlength)
 129             {
 130               sq.append(GapChar);
 131             }
 132           }
 133           else
 134           {
 135             sq.append(reference.substring(cursor, sbend));
 136           }
 137         }
 138         alcursor += range[i];
 139         cursor += range[i];
 140         end = cursor - 1;
 141         endpos = alcursor;
 142         break;
 143       default:
 144         throw new Error("Unknown SeqCigar operation '" + operation[i] + "'");
 145       }
 146     }
 147     if (++delcount > 0)
 148     {
 149       trunc_deletions = new int[delcount][];
 150       System.arraycopy(deletions, 0, trunc_deletions, 0, delcount);
 151     }
 152     deletions = null;
 153     return new Object[]
 154     { ((reference != null) ? sq.toString() : null), new int[]
 155     { start, startpos, end, endpos }, trunc_deletions };
 156   }
 157
 158   protected void compact_operations()
 159   {
 160     int i = 1;
 161     if (operation == null)
 162     {
 163       return;
 164     }
 165     char last = operation[0];
 166     while (i < length)
 167     {
 168       if (last == operation[i])
 169       {
 170         range[i - 1] += range[i];
 171         int r = length - i;
 172         if (r > 0)
 173         {
 174           System.arraycopy(range, i + 1, range, i, r);
 175           System.arraycopy(operation, i + 1, operation, i, r);
 176         }
 177         length--;
 178       }
 179       else
 180       {
 181         last = operation[i++];
 182       }
 183     }
 184   }
 185
 186   /**
 187    * turn a cigar string into a series of operation range pairs
 188    *
 189    * @param cigarString
 190    *          String
 191    * @return object[] {char[] operation, int[] range}
 192    * @throws java.lang.Exception
 193    *           for improperly formated cigar strings or ones with unknown
 194    *           operations
 195    */
 196   public static Object[] parseCigarString(String cigarString)
 197           throws Exception
 198   {
 199     int ops = 0;
 200     for (int i = 0, l = cigarString.length(); i < l; i++)
 201     {
 202       char c = cigarString.charAt(i);
 203       if (c == M || c == (M - _case_shift) || c == I
 204               || c == (I - _case_shift) || c == D || c == (D - _case_shift))
 205       {
 206         ops++;
 207       }
 208     }
 209     char[] operation = new char[ops];
 210     int[] range = new int[ops];
 211     int op = 0;
 212     int i = 0, l = cigarString.length();
 213     while (i < l)
 214     {
 215       char c;
 216       int j = i;
 217       do
 218       {
 219         c = cigarString.charAt(j++);
 220       } while (c >= '0' && c <= '9' && j < l);
 221       if (j >= l && c >= '0' && c <= '9')
 222       {
 223         throw new Exception("Unterminated cigar string.");
 224       }
 225       try
 226       {
 227         String rangeint = cigarString.substring(i, j - 1);
 228         range[op] = Integer.parseInt(rangeint);
 229         i = j;
 230       } catch (Exception e)
 231       {
 232         throw new Error("Implementation bug in parseCigarString");
 233       }
 234       if (c >= 'a' && c <= 'z')
 235       {
 236         c -= _case_shift;
 237       }
 238       if ((c == M || c == I || c == D))
 239       {
 240         operation[op++] = c;
 241       }
 242       else
 243       {
 244         throw new Exception("Unexpected operation '" + c
 245                 + "' in cigar string (position " + i + " in '"
 246                 + cigarString + "'");
 247       }
 248     }
 249     return new Object[]
 250     { operation, range };
 251   }
 252
 253   /**
 254    * add an operation to cigar string
 255    *
 256    * @param op
 257    *          char
 258    * @param range
 259    *          int
 260    */
 261   public void addOperation(char op, int range)
 262   {
 263     if (op >= 'a' && op <= 'z')
 264     {
 265       op -= _case_shift;
 266     }
 267     if (op != M && op != D && op != I)
 268     {
 269       throw new Error("Implementation error. Invalid operation string.");
 270     }
 271     if (range == 0)
 272     {
 273       return; // No Operation to add.
 274     }
 275     if (range < 0)
 276     {
 277       throw new Error(
 278               "Invalid range string (must be zero or positive number)");
 279     }
 280     int lngth = 0;
 281     if (operation == null)
 282     {
 283       this.operation = new char[_inc_length];
 284       this.range = new int[_inc_length];
 285     }
 286     if (length + 1 == operation.length)
 287     {
 288       char[] ops = this.operation;
 289       this.operation = new char[length + 1 + _inc_length];
 290       System.arraycopy(ops, 0, this.operation, 0, length);
 291       ops = null;
 292       int[] rng = this.range;
 293       this.range = new int[length + 1 + _inc_length];
 294       System.arraycopy(rng, 0, this.range, 0, length);
 295       rng = null;
 296     }
 297     if ((length > 0) && (operation[length - 1] == op))
 298     {
 299       length--; // modify existing operation.
 300     }
 301     else
 302     {
 303       this.range[length] = 0; // reset range
 304     }
 305     this.operation[length] = op;
 306     this.range[length++] += range;
 307   }
 308
 309   /**
 310    * semi-efficient insert an operation on the current cigar string set at
 311    * column pos (from 1) NOTE: Insertion operations simply extend width of cigar
 312    * result - affecting registration of alignment Deletion ops will shorten
 313    * length of result - and affect registration of alignment Match ops will also
 314    * affect length of result - affecting registration of alignment (ie
 315    * "10M".insert(4,I,3)->"4M3I3M") - (replace?) (ie
 316    * "10M".insert(4,D,3)->"4M3D3M") - (shortens alignment) (ie
 317    * "5I5M".insert(4,I,3)->"8I5M") - real insertion (ie
 318    * "5I5M".insert(4,D,3)->"4I2D3M") - shortens aligment - I's are removed, Ms
 319    * changed to Ds (ie "10M".insert(4,M,3)->"13M") - lengthens - Is changed to
 320    * M, Ds changed to M. (ie "5I5M".insert(4,M,3)->"4I8M") - effectively shifts
 321    * sequence left by 1 residue and extends it by 3 (
 322    * "10D5M".insert(-1,M,3)->"3M7D5M") ( "10D5M".insert(0,M,3)->"7D8M") (
 323    * "10D5M".insert(1,M,3)->"10D8M") ( "1M10D5M".insert(0,M,3)->"1M10D8M") (
 324    * "1M10D5M".insert(1,M,3)->"
 325    *
 326    * if pos is beyond width - I operations are added before the operation
 327    *
 328    * @param pos
 329    *          int -1, 0-length of visible region, or greater to append new ops
 330    *          (with insertions in between)
 331    * @param op
 332    *          char
 333    * @param range
 334    *          int public void addOperationAt(int pos, char op, int range) { int
 335    *          cursor = -1; // mark the position for the current operation being
 336    *          edited. int o = 0; boolean last_d = false; // previous op was a
 337    *          deletion. if (pos < -1) throw new
 338    *          Error("pos<-1 is not supported."); while (o<length) { if
 339    *          (operation[o] != D) { if ( (cursor + this.range[o]) < pos) {
 340    *          cursor += this.range[o]; o++; last_d=false; } else { break; } }
 341    *          else { last_d=true; o++; } } if (o==length) { // must insert more
 342    *          operations before pos if (pos-cursor>0) addInsertion(pos-cursor);
 343    *          // then just add the new operation. Regardless of what it is.
 344    *          addOperation(op, range); } else { int diff = pos - cursor;
 345    *
 346    *          int e_length = length-o; // new edit operation array length. //
 347    *          diff<0 - can only happen before first insertion or match. -
 348    *          affects op and all following // dif==0 - only when at first
 349    *          position of existing op - // diff>0 - must preserve some existing
 350    *          operations int[] e_range = new int[e_length];
 351    *          System.arraycopy(this.range, o, e_range, 0, e_length); char[] e_op
 352    *          = new char[e_length]; System.arraycopy(this.operation, o, e_op, 0,
 353    *          e_length); length = o; // can now use add_operation to extend
 354    *          list. int e_o=0; // current operation being edited. switch (op) {
 355    *          case M: switch (e_op[e_o]) { case M: if (last_d && diff <= 0) { //
 356    *          reduce D's, if possible if (range<=this.range[o-1]) { this.range[o
 357    *          - 1] -= range; } else { this.range[o-1]=0; } if
 358    *          (this.range[o-1]==0) o--; // lose this op. } e_range[e_o] +=
 359    *          range; // just add more matched residues break; case I: // change
 360    *          from insertion to match if (last_d && diff<=0) { // reduce D's, if
 361    *          possible if (range<=this.range[o-1]) { this.range[o - 1] -= range;
 362    *          } else { this.range[o-1]=0; } if (this.range[o-1]==0) o--; // lose
 363    *          this op. } e_range[e_o] break; default: throw new Inp }
 364    *
 365    *          break; case I: break; case D: } break; default: throw new
 366    *          Error("Implementation Error: Unknown operation in addOperation!");
 367    *          } // finally, add remaining ops. while (e_o<e_length) {
 368    *          addOperation(e_op[e_o], e_range[e_o]); e_o++; } } }
 369    */
 370   /**
 371    * Mark residues from start to end (inclusive) as deleted from the alignment,
 372    * and removes any insertions.
 373    *
 374    * @param start
 375    *          int
 376    * @param end
 377    *          int
 378    * @return deleted int - number of symbols marked as deleted
 379    */
 380   public int deleteRange(int start, int end)
 381   {
 382     int deleted = 0;
 383     if (length == 0)
 384     {
 385       // nothing to do here
 386       return deleted;
 387     }
 388     if (start < 0 || start > end)
 389     {
 390       throw new Error(
 391               "Implementation Error: deleteRange out of bounds: start must be non-negative and less than end.");
 392     }
 393     // find beginning
 394     int cursor = 0; // mark the position for the current operation being edited.
 395     int rlength = 1 + end - start; // number of positions to delete
 396     int oldlen = length;
 397     int o = 0;
 398     boolean editing = false;
 399     char[] oldops = operation;
 400     int[] oldrange = range;
 401     length = 0;
 402     operation = null;
 403     range = null;
 404     compact_operations();
 405     while (o < oldlen && cursor <= end && rlength > 0)
 406     {
 407       if (oldops[o] == D)
 408       {
 409         // absorbed into new deleted region.
 410         addDeleted(oldrange[o++]);
 411         continue;
 412       }
 413
 414       int remain = oldrange[o]; // number of op characters left to edit
 415       if (!editing)
 416       {
 417         if ((cursor + remain) <= start)
 418         {
 419           addOperation(oldops[o], oldrange[o]);
 420           cursor += oldrange[o++];
 421           continue; // next operation
 422         }
 423         editing = true;
 424         // add operations before hidden region
 425         if (start - cursor > 0)
 426         {
 427           addOperation(oldops[o], start - cursor);
 428           remain -= start - cursor;
 429         }
 430       }
 431       // start inserting new ops
 432       if (o < oldlen && editing && rlength > 0 && remain > 0)
 433       {
 434         switch (oldops[o])
 435         {
 436         case M:
 437           if (rlength > remain)
 438           {
 439             addDeleted(remain);
 440             deleted += remain;
 441           }
 442           else
 443           {
 444             deleted += rlength;
 445             addDeleted(rlength);
 446             if (remain - rlength > 0)
 447             {
 448               this.addOperation(M, remain - rlength); // add remaining back.
 449             }
 450             rlength = 0;
 451             remain = 0;
 452           }
 453           break;
 454         case I:
 455           if (remain - rlength > 0)
 456           {
 457             // only remove some gaps
 458             addInsertion(remain - rlength);
 459             rlength = 0;
 460           }
 461           break;
 462         case D:
 463           throw new Error("Implementation error."); // do nothing;
 464         default:
 465           throw new Error("Implementation Error! Unknown operation '"
 466                   + oldops[o] + "'");
 467         }
 468         rlength -= remain;
 469         remain = oldrange[++o]; // number of op characters left to edit
 470       }
 471     }
 472     // add remaining
 473     while (o < oldlen)
 474     {
 475       addOperation(oldops[o], oldrange[o++]);
 476     }
 477     // if (cursor<(start+1)) {
 478     // ran out of ops - nothing to do here ?
 479     // addInsertion(start-cursor);
 480     // }
 481     return deleted;
 482   }
 483
 484   /**
 485    * Deleted regions mean that there will be discontinuous sequence numbering in
 486    * the sequence returned by getSeq(char).
 487    *
 488    * @return true if there deletions
 489    */
 490   public boolean hasDeletedRegions()
 491   {
 492     for (int i = 0; i < length; i++)
 493     {
 494       if (operation[i] == D)
 495       {
 496         return true;
 497       }
 498     }
 499     return false;
 500   }
 501
 502   /**
 503    * enumerate the ranges on seq that are marked as deleted in this cigar
 504    *
 505    * @return int[] { vis_start, sym_start, length }
 506    */
 507   public int[] getDeletedRegions()
 508   {
 509     if (length == 0)
 510     {
 511       return null;
 512     }
 513     Vector dr = new Vector();
 514     int cursor = 0, vcursor = 0;
 515     for (int i = 0; i < length; i++)
 516     {
 517       switch (operation[i])
 518       {
 519       case M:
 520         cursor += range[i];
 521       case I:
 522         vcursor += range[i];
 523         break;
 524       case D:
 525         dr.addElement(new int[]
 526         { vcursor, cursor, range[i] });
 527         cursor += range[i];
 528       }
 529     }
 530     if (dr.size() == 0)
 531     {
 532       return null;
 533     }
 534     int[] delregions = new int[dr.size() * 3];
 535     for (int i = 0, l = dr.size(); i < l; i++)
 536     {
 537       int[] reg = (int[]) dr.elementAt(i);
 538       delregions[i * 3] = reg[0];
 539       delregions[i * 3 + 1] = reg[1];
 540       delregions[i * 3 + 2] = reg[2];
 541     }
 542     return delregions;
 543   }
 544
 545   /**
 546    * sum of ranges in cigar string
 547    *
 548    * @return int number of residues hidden, matched, or gaps inserted into
 549    *         sequence
 550    */
 551   public int getFullWidth()
 552   {
 553     int w = 0;
 554     if (range != null)
 555     {
 556       for (int i = 0; i < length; i++)
 557       {
 558         w += range[i];
 559       }
 560     }
 561     return w;
 562   }
 563
 564   /**
 565    * Visible length of aligned sequence
 566    *
 567    * @return int length of including gaps and less hidden regions
 568    */
 569   public int getWidth()
 570   {
 571     int w = 0;
 572     if (range != null)
 573     {
 574       for (int i = 0; i < length; i++)
 575       {
 576         if (operation[i] == M || operation[i] == I)
 577         {
 578           w += range[i];
 579         }
 580       }
 581     }
 582     return w;
 583   }
 584
 585   /**
 586    * mark a range of inserted residues
 587    *
 588    * @param range
 589    *          int
 590    */
 591   public void addInsertion(int range)
 592   {
 593     this.addOperation(I, range);
 594   }
 595
 596   /**
 597    * mark the next range residues as hidden (not aligned) or deleted
 598    *
 599    * @param range
 600    *          int
 601    */
 602   public void addDeleted(int range)
 603   {
 604     this.addOperation(D, range);
 605   }
 606
 607   /**
 608    * Modifies operation list to delete columns from start to end (inclusive)
 609    * editing will remove insertion operations, and convert matches to deletions
 610    *
 611    * @param start
 612    *          alignment column
 613    * @param end
 614    *          alignment column
 615    * @return boolean true if residues were marked as deleted. public boolean
 616    *         deleteRange(int start, int end) { boolean deleted = false; int op =
 617    *         0, prevop = -1, firstm = -1, lastm = -1, postop = -1; int width =
 618    *         0; // zero'th column if (length > 0) { // find operation bracketing
 619    *         start of the range do { if (operation[op] != D) { width +=
 620    *         range[prevop = op]; } op++; } while (op < length && width < start);
 621    *         } if (width < start) { // run off end - add more operations up to
 622    *         deletion. addInsertion(start - width); } else { // edit existing
 623    *         operations. op = prevop; width -= range[prevop]; int[] oldrange =
 624    *         range; char[] oldops = operation; range = new int[oldrange.length];
 625    *         operation = new char[oldops.length]; if (op < length) { do { if
 626    *         (operation[op] != D) { width += range[postop = op]; } op++; } while
 627    *         (op < length && width <= end); } } if (deleted == true) {
 628    *         addDeleted(end - start + 1); } return deleted; }
 629    */
 630   /**
 631    * Return an ENSEMBL style cigar string where D may indicates excluded parts
 632    * of seq
 633    *
 634    * @return String of form ([0-9]+[IMD])+
 635    */
 636   public String getCigarstring()
 637   {
 638     StringBuffer cigarString = new StringBuffer();
 639     for (int i = 0; i < length; i++)
 640     {
 641       cigarString.append("" + range[i]);
 642       cigarString.append(operation[i]);
 643     }
 644     return cigarString.toString();
 645   }
 646 }