src/jalview/datamodel/CigarBase.java

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