JAL-1705 include any unmapped protein start 'X' when aligning to dna
[jalview.git] / src / jalview / datamodel / Mapping.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.MapList;
24
25 import java.util.Iterator;
26 import java.util.NoSuchElementException;
27 import java.util.Vector;
28
29 public class Mapping
30 {
31   /**
32    * An iterator that serves the aligned codon positions (with their protein
33    * products).
34    * 
35    * @author gmcarstairs
36    *
37    */
38   public class AlignedCodonIterator implements Iterator<AlignedCodon>
39   {
40     /*
41      * The gap character used in the aligned sequence
42      */
43     private final char gap;
44
45     /*
46      * The characters of the aligned sequence e.g. "-cGT-ACgTG-"
47      */
48     private final char[] alignedSeq;
49
50     /*
51      * the sequence start residue
52      */
53     private int start;
54
55     /*
56      * Next position (base 0) in the aligned sequence
57      */
58     private int alignedColumn = 0;
59
60     /*
61      * Count of bases up to and including alignedColumn position
62      */
63     private int alignedBases = 0;
64
65     /*
66      * [start, end] from ranges (base 1)
67      */
68     private Iterator<int[]> fromRanges;
69
70     /*
71      * [start, end] to ranges (base 1)
72      */
73     private Iterator<int[]> toRanges;
74
75     /*
76      * The current [start, end] (base 1) from range
77      */
78     private int[] currentFromRange = null;
79
80     /*
81      * The current [start, end] (base 1) to range
82      */
83     private int[] currentToRange = null;
84
85     /*
86      * The next 'from' position (base 1) to process
87      */
88     private int fromPosition = 0;
89
90     /*
91      * The next 'to' position (base 1) to process
92      */
93     private int toPosition = 0;
94
95     /**
96      * Constructor
97      * 
98      * @param seq
99      *          the aligned sequence
100      * @param gapChar
101      */
102     public AlignedCodonIterator(SequenceI seq, char gapChar)
103     {
104       this.alignedSeq = seq.getSequence();
105       this.start = seq.getStart();
106       this.gap = gapChar;
107       fromRanges = map.getFromRanges().iterator();
108       toRanges = map.getToRanges().iterator();
109       if (fromRanges.hasNext())
110       {
111         currentFromRange = fromRanges.next();
112         fromPosition = currentFromRange[0];
113       }
114       if (toRanges.hasNext())
115       {
116         currentToRange = toRanges.next();
117         toPosition = currentToRange[0];
118       }
119     }
120
121     /**
122      * Returns true unless we have already traversed the whole mapping.
123      */
124     @Override
125     public boolean hasNext()
126     {
127       if (fromRanges.hasNext())
128       {
129         return true;
130       }
131       if (currentFromRange == null || fromPosition >= currentFromRange[1])
132       {
133         return false;
134       }
135       return true;
136     }
137
138     /**
139      * Returns the next codon's aligned positions, and translated value.
140      * 
141      * @throws NoSuchElementException
142      *           if hasNext() would have returned false
143      * @throws IncompleteCodonException
144      *           if not enough mapped bases are left to make up a codon
145      */
146     @Override
147     public AlignedCodon next() throws IncompleteCodonException
148     {
149       if (!hasNext())
150       {
151         throw new NoSuchElementException();
152       }
153
154       int[] codon = getNextCodon();
155       int[] alignedCodon = getAlignedCodon(codon);
156
157       String peptide = getPeptide();
158       int peptideCol = toPosition - 1 - Mapping.this.to.getStart();
159       return new AlignedCodon(alignedCodon[0], alignedCodon[1],
160               alignedCodon[2], peptide, peptideCol);
161     }
162
163     /**
164      * Retrieve the translation as the 'mapped to' position in the mapped to
165      * sequence.
166      * 
167      * @return
168      * @throws NoSuchElementException
169      *           if the 'toRange' is exhausted (nothing to map to)
170      */
171     private String getPeptide()
172     {
173       // TODO should ideally handle toRatio other than 1 as well...
174       // i.e. code like getNextCodon()
175       if (toPosition <= currentToRange[1])
176       {
177         SequenceI seq = Mapping.this.to;
178         char pep = seq.getSequence()[toPosition - seq.getStart()];
179         toPosition++;
180         return String.valueOf(pep);
181       }
182       if (!toRanges.hasNext())
183       {
184         throw new NoSuchElementException("Ran out of peptide at position "
185                 + toPosition);
186       }
187       currentToRange = toRanges.next();
188       toPosition = currentToRange[0];
189       return getPeptide();
190     }
191
192     /**
193      * Get the (base 1) dataset positions for the next codon in the mapping.
194      * 
195      * @throws IncompleteCodonException
196      *           if less than 3 remaining bases are mapped
197      */
198     private int[] getNextCodon()
199     {
200       int[] codon = new int[3];
201       int codonbase = 0;
202
203       while (codonbase < 3)
204       {
205         if (fromPosition <= currentFromRange[1])
206         {
207           /*
208            * Add next position from the current start-end range
209            */
210           codon[codonbase++] = fromPosition++;
211         }
212         else
213         {
214           /*
215            * Move to the next range - if there is one
216            */
217           if (!fromRanges.hasNext())
218           {
219             throw new IncompleteCodonException();
220           }
221           currentFromRange = fromRanges.next();
222           fromPosition = currentFromRange[0];
223         }
224       }
225       return codon;
226     }
227
228     /**
229      * Get the aligned column positions (base 0) for the given sequence
230      * positions (base 1), by counting ungapped characters in the aligned
231      * sequence.
232      * 
233      * @param codon
234      * @return
235      */
236     private int[] getAlignedCodon(int[] codon)
237     {
238       int[] aligned = new int[codon.length];
239       for (int i = 0; i < codon.length; i++)
240       {
241         aligned[i] = getAlignedColumn(codon[i]);
242       }
243       return aligned;
244     }
245
246     /**
247      * Get the aligned column position (base 0) for the given sequence position
248      * (base 1).
249      * 
250      * @param sequencePos
251      * @return
252      */
253     private int getAlignedColumn(int sequencePos)
254     {
255       /*
256        * allow for offset e.g. treat pos 8 as 2 if sequence starts at 7
257        */
258       int truePos = sequencePos - (start - 1);
259       while (alignedBases < truePos && alignedColumn < alignedSeq.length)
260       {
261         if (alignedSeq[alignedColumn++] != gap)
262         {
263           alignedBases++;
264         }
265       }
266       return alignedColumn - 1;
267     }
268
269     @Override
270     public void remove()
271     {
272       // ignore
273     }
274
275   }
276
277   /**
278    * Contains the start-end pairs mapping from the associated sequence to the
279    * sequence in the database coordinate system. It also takes care of step
280    * difference between coordinate systems.
281    */
282   MapList map = null;
283
284   /**
285    * The sequence that map maps the associated sequence to (if any).
286    */
287   SequenceI to = null;
288
289   public Mapping(MapList map)
290   {
291     super();
292     this.map = map;
293   }
294
295   public Mapping(SequenceI to, MapList map)
296   {
297     this(map);
298     this.to = to;
299   }
300
301   /**
302    * create a new mapping from
303    * 
304    * @param to
305    *          the sequence being mapped
306    * @param exon
307    *          int[] {start,end,start,end} series on associated sequence
308    * @param is
309    *          int[] {start,end,...} ranges on the reference frame being mapped
310    *          to
311    * @param i
312    *          step size on associated sequence
313    * @param j
314    *          step size on mapped frame
315    */
316   public Mapping(SequenceI to, int[] exon, int[] is, int i, int j)
317   {
318     this(to, new MapList(exon, is, i, j));
319   }
320
321   /**
322    * create a duplicate (and independent) mapping object with the same reference
323    * to any SequenceI being mapped to.
324    * 
325    * @param map2
326    */
327   public Mapping(Mapping map2)
328   {
329     if (map2 != this && map2 != null)
330     {
331       if (map2.map != null)
332       {
333         map = new MapList(map2.map);
334       }
335       to = map2.to;
336     }
337   }
338
339   /**
340    * @return the map
341    */
342   public MapList getMap()
343   {
344     return map;
345   }
346
347   /**
348    * @param map
349    *          the map to set
350    */
351   public void setMap(MapList map)
352   {
353     this.map = map;
354   }
355
356   /**
357    * Equals that compares both the to references and MapList mappings.
358    * 
359    * @param other
360    * @return
361    */
362   @Override
363   public boolean equals(Object o)
364   {
365     // TODO should override Object.hashCode() to ensure that equal objects have
366     // equal hashcodes
367     if (o == null || !(o instanceof Mapping))
368     {
369       return false;
370     }
371     Mapping other = (Mapping) o;
372     if (other == this)
373     {
374       return true;
375     }
376     if (other.to != to)
377     {
378       return false;
379     }
380     if ((map != null && other.map == null)
381             || (map == null && other.map != null))
382     {
383       return false;
384     }
385     if ((map == null && other.map == null) || map.equals(other.map))
386     {
387       return true;
388     }
389     return false;
390   }
391
392   /**
393    * get the 'initial' position in the associated sequence for a position in the
394    * mapped reference frame
395    * 
396    * @param mpos
397    * @return
398    */
399   public int getPosition(int mpos)
400   {
401     if (map != null)
402     {
403       int[] mp = map.shiftTo(mpos);
404       if (mp != null)
405       {
406         return mp[0];
407       }
408     }
409     return mpos;
410   }
411
412   /**
413    * gets boundary in direction of mapping
414    * 
415    * @param position
416    *          in mapped reference frame
417    * @return int{start, end} positions in associated sequence (in direction of
418    *         mapped word)
419    */
420   public int[] getWord(int mpos)
421   {
422     if (map != null)
423     {
424       return map.getToWord(mpos);
425     }
426     return null;
427   }
428
429   /**
430    * width of mapped unit in associated sequence
431    * 
432    */
433   public int getWidth()
434   {
435     if (map != null)
436     {
437       return map.getFromRatio();
438     }
439     return 1;
440   }
441
442   /**
443    * width of unit in mapped reference frame
444    * 
445    * @return
446    */
447   public int getMappedWidth()
448   {
449     if (map != null)
450     {
451       return map.getToRatio();
452     }
453     return 1;
454   }
455
456   /**
457    * get mapped position in the associated reference frame for position pos in
458    * the associated sequence.
459    * 
460    * @param pos
461    * @return
462    */
463   public int getMappedPosition(int pos)
464   {
465     if (map != null)
466     {
467       int[] mp = map.shiftFrom(pos);
468       if (mp != null)
469       {
470         return mp[0];
471       }
472     }
473     return pos;
474   }
475
476   public int[] getMappedWord(int pos)
477   {
478     if (map != null)
479     {
480       int[] mp = map.shiftFrom(pos);
481       if (mp != null)
482       {
483         return new int[] { mp[0], mp[0] + mp[2] * (map.getToRatio() - 1) };
484       }
485     }
486     return null;
487   }
488
489   /**
490    * locates the region of feature f in the associated sequence's reference
491    * frame
492    * 
493    * @param f
494    * @return one or more features corresponding to f
495    */
496   public SequenceFeature[] locateFeature(SequenceFeature f)
497   {
498     if (true)
499     { // f.getBegin()!=f.getEnd()) {
500       if (map != null)
501       {
502         int[] frange = map.locateInFrom(f.getBegin(), f.getEnd());
503         if (frange == null)
504         {
505           // JBPNote - this isprobably not the right thing to doJBPHack
506           return null;
507         }
508         SequenceFeature[] vf = new SequenceFeature[frange.length / 2];
509         for (int i = 0, v = 0; i < frange.length; i += 2, v++)
510         {
511           vf[v] = new SequenceFeature(f);
512           vf[v].setBegin(frange[i]);
513           vf[v].setEnd(frange[i + 1]);
514           if (frange.length > 2)
515           {
516             vf[v].setDescription(f.getDescription() + "\nPart " + (v + 1));
517           }
518         }
519         return vf;
520       }
521     }
522     if (false) // else
523     {
524       int[] word = getWord(f.getBegin());
525       if (word[0] < word[1])
526       {
527         f.setBegin(word[0]);
528       }
529       else
530       {
531         f.setBegin(word[1]);
532       }
533       word = getWord(f.getEnd());
534       if (word[0] > word[1])
535       {
536         f.setEnd(word[0]);
537       }
538       else
539       {
540         f.setEnd(word[1]);
541       }
542     }
543     // give up and just return the feature.
544     return new SequenceFeature[] { f };
545   }
546
547   /**
548    * return a series of contigs on the associated sequence corresponding to the
549    * from,to interval on the mapped reference frame
550    * 
551    * @param from
552    * @param to
553    * @return int[] { from_i, to_i for i=1 to n contiguous regions in the
554    *         associated sequence}
555    */
556   public int[] locateRange(int from, int to)
557   {
558     if (map != null)
559     {
560       if (from <= to)
561       {
562         from = (map.getToLowest() < from) ? from : map.getToLowest();
563         to = (map.getToHighest() > to) ? to : map.getToHighest();
564         if (from > to)
565         {
566           return null;
567         }
568       }
569       else
570       {
571         from = (map.getToHighest() > from) ? from : map.getToHighest();
572         to = (map.getToLowest() < to) ? to : map.getToLowest();
573         if (from < to)
574         {
575           return null;
576         }
577       }
578       return map.locateInFrom(from, to);
579     }
580     return new int[] { from, to };
581   }
582
583   /**
584    * return a series of mapped contigs mapped from a range on the associated
585    * sequence
586    * 
587    * @param from
588    * @param to
589    * @return
590    */
591   public int[] locateMappedRange(int from, int to)
592   {
593     if (map != null)
594     {
595
596       if (from <= to)
597       {
598         from = (map.getFromLowest() < from) ? from : map.getFromLowest();
599         to = (map.getFromHighest() > to) ? to : map.getFromHighest();
600         if (from > to)
601         {
602           return null;
603         }
604       }
605       else
606       {
607         from = (map.getFromHighest() > from) ? from : map.getFromHighest();
608         to = (map.getFromLowest() < to) ? to : map.getFromLowest();
609         if (from < to)
610         {
611           return null;
612         }
613       }
614       return map.locateInTo(from, to);
615     }
616     return new int[] { from, to };
617   }
618
619   /**
620    * return a new mapping object with a maplist modifed to only map the visible
621    * regions defined by viscontigs.
622    * 
623    * @param viscontigs
624    * @return
625    */
626   public Mapping intersectVisContigs(int[] viscontigs)
627   {
628     Mapping copy = new Mapping(this);
629     if (map != null)
630     {
631       int vpos = 0;
632       int apos = 0;
633       Vector toRange = new Vector();
634       Vector fromRange = new Vector();
635       for (int vc = 0; vc < viscontigs.length; vc += 2)
636       {
637         // find a mapped range in this visible region
638         int[] mpr = locateMappedRange(1 + viscontigs[vc],
639                 viscontigs[vc + 1] - 1);
640         if (mpr != null)
641         {
642           for (int m = 0; m < mpr.length; m += 2)
643           {
644             toRange.addElement(new int[] { mpr[m], mpr[m + 1] });
645             int[] xpos = locateRange(mpr[m], mpr[m + 1]);
646             for (int x = 0; x < xpos.length; x += 2)
647             {
648               fromRange.addElement(new int[] { xpos[x], xpos[x + 1] });
649             }
650           }
651         }
652       }
653       int[] from = new int[fromRange.size() * 2];
654       int[] to = new int[toRange.size() * 2];
655       int[] r;
656       for (int f = 0, fSize = fromRange.size(); f < fSize; f++)
657       {
658         r = (int[]) fromRange.elementAt(f);
659         from[f * 2] = r[0];
660         from[f * 2 + 1] = r[1];
661       }
662       for (int f = 0, fSize = toRange.size(); f < fSize; f++)
663       {
664         r = (int[]) toRange.elementAt(f);
665         to[f * 2] = r[0];
666         to[f * 2 + 1] = r[1];
667       }
668       copy.setMap(new MapList(from, to, map.getFromRatio(), map
669               .getToRatio()));
670     }
671     return copy;
672   }
673
674   /**
675    * get the sequence being mapped to - if any
676    * 
677    * @return null or a dataset sequence
678    */
679   public SequenceI getTo()
680   {
681     return to;
682   }
683
684   /**
685    * set the dataset sequence being mapped to if any
686    * 
687    * @param tto
688    */
689   public void setTo(SequenceI tto)
690   {
691     to = tto;
692   }
693
694   /*
695    * (non-Javadoc)
696    * 
697    * @see java.lang.Object#finalize()
698    */
699   @Override
700   protected void finalize() throws Throwable
701   {
702     map = null;
703     to = null;
704     super.finalize();
705   }
706
707   /**
708    * Returns an iterator which can serve up the aligned codon column positions
709    * and their corresponding peptide products
710    * 
711    * @param seq
712    *          an aligned (i.e. possibly gapped) sequence
713    * @param gapChar
714    * @return
715    */
716   public Iterator<AlignedCodon> getCodonIterator(SequenceI seq, char gapChar)
717   {
718     return new AlignedCodonIterator(seq, gapChar);
719   }
720
721   /**
722    * Readable representation for debugging only, not guaranteed not to change
723    */
724   @Override
725   public String toString()
726   {
727     return String.format("%s %s", this.map.toString(), this.to == null ? ""
728             : this.to.getName());
729   }
730
731 }