src/jalview/datamodel/CigarBase.java

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