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