src/jalview/datamodel/CigarBase.java

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