src/jalview/analysis/Dna.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.analysis;
  22
  23 import jalview.api.AlignViewportI;
  24 import jalview.bin.Cache;
  25 import jalview.datamodel.AlignedCodonFrame;
  26 import jalview.datamodel.Alignment;
  27 import jalview.datamodel.AlignmentAnnotation;
  28 import jalview.datamodel.AlignmentI;
  29 import jalview.datamodel.Annotation;
  30 import jalview.datamodel.DBRefEntry;
  31 import jalview.datamodel.DBRefSource;
  32 import jalview.datamodel.FeatureProperties;
  33 import jalview.datamodel.GraphLine;
  34 import jalview.datamodel.Mapping;
  35 import jalview.datamodel.Sequence;
  36 import jalview.datamodel.SequenceFeature;
  37 import jalview.datamodel.SequenceI;
  38 import jalview.schemes.ResidueProperties;
  39 import jalview.util.Comparison;
  40 import jalview.util.DBRefUtils;
  41 import jalview.util.MapList;
  42 import jalview.util.ShiftList;
  43
  44 import java.util.ArrayList;
  45 import java.util.Arrays;
  46 import java.util.Hashtable;
  47 import java.util.List;
  48
  49 public class Dna
  50 {
  51   private static final String STOP_X = "X";
  52
  53   final private List<SequenceI> selection;
  54
  55   final private String[] seqstring;
  56
  57   final private int[] contigs;
  58
  59   final private char gapChar;
  60
  61   final private AlignmentAnnotation[] annotations;
  62
  63   final private int dnaWidth;
  64
  65   final private Alignment dataset;
  66
  67   private int aaWidth = 0;
  68
  69   /**
  70    * Constructor given a viewport and the visible contigs.
  71    *
  72    * @param viewport
  73    * @param visibleContigs
  74    */
  75   public Dna(AlignViewportI viewport, int[] visibleContigs)
  76   {
  77     this.selection = Arrays.asList(viewport.getSequenceSelection());
  78     this.seqstring = viewport.getViewAsString(true);
  79     this.contigs = visibleContigs;
  80     this.gapChar = viewport.getGapCharacter();
  81     this.annotations = viewport.getAlignment().getAlignmentAnnotation();
  82     this.dnaWidth = viewport.getAlignment().getWidth();
  83     this.dataset = viewport.getAlignment().getDataset();
  84   }
  85
  86   /**
  87    * Test whether codon positions cdp1 should align before, with, or after cdp2.
  88    * Returns zero if all positions match (or either argument is null). Returns
  89    * -1 if any position in the first codon precedes the corresponding position
  90    * in the second codon. Else returns +1 (some position in the second codon
  91    * precedes the corresponding position in the first).
  92    *
  93    * Note this is not necessarily symmetric, for example:
  94    * <ul>
  95    * <li>compareCodonPos([2,5,6], [3,4,5]) returns -1</li>
  96    * <li>compareCodonPos([3,4,5], [2,5,6]) also returns -1</li>
  97    * </ul>
  98    *
  99    * @param cdp1
 100    * @param cdp2
 101    * @return
 102    */
 103   public static int compareCodonPos(int[] cdp1, int[] cdp2)
 104   {
 105     if (cdp1 == null
 106             || cdp2 == null
 107             || (cdp1[0] == cdp2[0] && cdp1[1] == cdp2[1] && cdp1[2] == cdp2[2]))
 108     {
 109       return 0;
 110     }
 111     if (cdp1[0] < cdp2[0] || cdp1[1] < cdp2[1] || cdp1[2] < cdp2[2])
 112     {
 113       // one base in cdp1 precedes the corresponding base in the other codon
 114       return -1;
 115     }
 116     // one base in cdp1 appears after the corresponding base in the other codon.
 117     return 1;
 118   }
 119
 120   /**
 121    *
 122    * @return
 123    */
 124   public AlignmentI translateCdna()
 125   {
 126     AlignedCodonFrame acf = new AlignedCodonFrame();
 127
 128     /*
 129      * This array will be built up so that position i holds the codon positions
 130      * e.g. [7, 9, 10] that match column i (base 0) in the aligned translation.
 131      * Note this implies a contract that if two codons do not align exactly,
 132      * their translated products must occupy different column positions.
 133      */
 134     int[][] alignedCodons = new int[dnaWidth][];
 135
 136     int s;
 137     int sSize = selection.size();
 138     List<SequenceI> pepseqs = new ArrayList<SequenceI>();
 139     for (s = 0; s < sSize; s++)
 140     {
 141       SequenceI newseq = translateCodingRegion(selection.get(s),
 142               seqstring[s], acf, alignedCodons, pepseqs);
 143
 144       if (newseq != null)
 145       {
 146         pepseqs.add(newseq);
 147         SequenceI ds = newseq;
 148         if (dataset != null)
 149         {
 150           while (ds.getDatasetSequence() != null)
 151           {
 152             ds = ds.getDatasetSequence();
 153           }
 154           dataset.addSequence(ds);
 155         }
 156       }
 157     }
 158
 159     SequenceI[] newseqs = pepseqs.toArray(new SequenceI[pepseqs.size()]);
 160     AlignmentI al = new Alignment(newseqs);
 161     // ensure we look aligned.
 162     al.padGaps();
 163     // link the protein translation to the DNA dataset
 164     al.setDataset(dataset);
 165     translateAlignedAnnotations(annotations, al, acf, alignedCodons);
 166     al.addCodonFrame(acf);
 167     return al;
 168   }
 169
 170   /**
 171    * fake the collection of DbRefs with associated exon mappings to identify if
 172    * a translation would generate distinct product in the currently selected
 173    * region.
 174    *
 175    * @param selection
 176    * @param viscontigs
 177    * @return
 178    */
 179   public static boolean canTranslate(SequenceI[] selection,
 180           int viscontigs[])
 181   {
 182     for (int gd = 0; gd < selection.length; gd++)
 183     {
 184       SequenceI dna = selection[gd];
 185       DBRefEntry[] dnarefs = DBRefUtils
 186               .selectRefs(dna.getDBRef(),
 187                       jalview.datamodel.DBRefSource.DNACODINGDBS);
 188       if (dnarefs != null)
 189       {
 190         // intersect with pep
 191         List<DBRefEntry> mappedrefs = new ArrayList<DBRefEntry>();
 192         DBRefEntry[] refs = dna.getDBRef();
 193         for (int d = 0; d < refs.length; d++)
 194         {
 195           if (refs[d].getMap() != null && refs[d].getMap().getMap() != null
 196                   && refs[d].getMap().getMap().getFromRatio() == 3
 197                   && refs[d].getMap().getMap().getToRatio() == 1)
 198           {
 199             mappedrefs.add(refs[d]); // add translated protein maps
 200           }
 201         }
 202         dnarefs = mappedrefs.toArray(new DBRefEntry[mappedrefs.size()]);
 203         for (int d = 0; d < dnarefs.length; d++)
 204         {
 205           Mapping mp = dnarefs[d].getMap();
 206           if (mp != null)
 207           {
 208             for (int vc = 0; vc < viscontigs.length; vc += 2)
 209             {
 210               int[] mpr = mp.locateMappedRange(viscontigs[vc],
 211                       viscontigs[vc + 1]);
 212               if (mpr != null)
 213               {
 214                 return true;
 215               }
 216             }
 217           }
 218         }
 219       }
 220     }
 221     return false;
 222   }
 223
 224   /**
 225    * Translate na alignment annotations onto translated amino acid alignment al
 226    * using codon mapping codons
 227    *
 228    * @param annotations
 229    * @param al
 230    * @param acf
 231    */
 232   protected static void translateAlignedAnnotations(
 233           AlignmentAnnotation[] annotations, AlignmentI al,
 234           AlignedCodonFrame acf, int[][] codons)
 235   {
 236     // Can only do this for columns with consecutive codons, or where
 237     // annotation is sequence associated.
 238
 239     if (annotations != null)
 240     {
 241       for (AlignmentAnnotation annotation : annotations)
 242       {
 243         /*
 244          * Skip hidden or autogenerated annotation. Also (for now), RNA
 245          * secondary structure annotation. If we want to show this against
 246          * protein we need a smarter way to 'translate' without generating
 247          * invalid (unbalanced) structure annotation.
 248          */
 249         if (annotation.autoCalculated || !annotation.visible
 250                 || annotation.isRNA())
 251         {
 252           continue;
 253         }
 254
 255         int aSize = acf.getaaWidth(); // aa alignment width.
 256         Annotation[] anots = (annotation.annotations == null) ? null
 257                 : new Annotation[aSize];
 258         if (anots != null)
 259         {
 260           for (int a = 0; a < aSize; a++)
 261           {
 262             // process through codon map.
 263             if (a < codons.length && codons[a] != null
 264                     && codons[a][0] == (codons[a][2] - 2))
 265             {
 266               anots[a] = getCodonAnnotation(codons[a],
 267                       annotation.annotations);
 268             }
 269           }
 270         }
 271
 272         AlignmentAnnotation aa = new AlignmentAnnotation(annotation.label,
 273                 annotation.description, anots);
 274         aa.graph = annotation.graph;
 275         aa.graphGroup = annotation.graphGroup;
 276         aa.graphHeight = annotation.graphHeight;
 277         if (annotation.getThreshold() != null)
 278         {
 279           aa.setThreshold(new GraphLine(annotation
 280                   .getThreshold()));
 281         }
 282         if (annotation.hasScore)
 283         {
 284           aa.setScore(annotation.getScore());
 285         }
 286
 287         final SequenceI seqRef = annotation.sequenceRef;
 288         if (seqRef != null)
 289         {
 290           SequenceI aaSeq = acf.getAaForDnaSeq(seqRef);
 291           if (aaSeq != null)
 292           {
 293             // aa.compactAnnotationArray(); // throw away alignment annotation
 294             // positioning
 295             aa.setSequenceRef(aaSeq);
 296             // rebuild mapping
 297             aa.createSequenceMapping(aaSeq, aaSeq.getStart(), true);
 298             aa.adjustForAlignment();
 299             aaSeq.addAlignmentAnnotation(aa);
 300           }
 301         }
 302         al.addAnnotation(aa);
 303       }
 304     }
 305   }
 306
 307   private static Annotation getCodonAnnotation(int[] is,
 308           Annotation[] annotations)
 309   {
 310     // Have a look at all the codon positions for annotation and put the first
 311     // one found into the translated annotation pos.
 312     int contrib = 0;
 313     Annotation annot = null;
 314     for (int p = 0; p < 3; p++)
 315     {
 316       if (annotations[is[p]] != null)
 317       {
 318         if (annot == null)
 319         {
 320           annot = new Annotation(annotations[is[p]]);
 321           contrib = 1;
 322         }
 323         else
 324         {
 325           // merge with last
 326           Annotation cpy = new Annotation(annotations[is[p]]);
 327           if (annot.colour == null)
 328           {
 329             annot.colour = cpy.colour;
 330           }
 331           if (annot.description == null || annot.description.length() == 0)
 332           {
 333             annot.description = cpy.description;
 334           }
 335           if (annot.displayCharacter == null)
 336           {
 337             annot.displayCharacter = cpy.displayCharacter;
 338           }
 339           if (annot.secondaryStructure == 0)
 340           {
 341             annot.secondaryStructure = cpy.secondaryStructure;
 342           }
 343           annot.value += cpy.value;
 344           contrib++;
 345         }
 346       }
 347     }
 348     if (contrib > 1)
 349     {
 350       annot.value /= contrib;
 351     }
 352     return annot;
 353   }
 354
 355   /**
 356    * Translate a na sequence
 357    *
 358    * @param selection
 359    *          sequence displayed under viscontigs visible columns
 360    * @param seqstring
 361    *          ORF read in some global alignment reference frame
 362    * @param viscontigs
 363    *          mapping from global reference frame to visible seqstring ORF read
 364    * @param acf
 365    *          Definition of global ORF alignment reference frame
 366    * @param alignedCodons
 367    * @param proteinSeqs
 368    * @param gapCharacter
 369    * @param starForStop
 370    *          when true stop codons will translate as '*', otherwise as 'X'
 371    * @return sequence ready to be added to alignment.
 372    */
 373   protected SequenceI translateCodingRegion(SequenceI selection,
 374           String seqstring, AlignedCodonFrame acf,
 375           int[][] alignedCodons, List<SequenceI> proteinSeqs)
 376   {
 377     List<int[]> skip = new ArrayList<int[]>();
 378     int skipint[] = null;
 379     ShiftList vismapping = new ShiftList(); // map from viscontigs to seqstring
 380     // intervals
 381     int vc;
 382     int[] scontigs = new int[contigs.length];
 383     int npos = 0;
 384     for (vc = 0; vc < contigs.length; vc += 2)
 385     {
 386       if (vc == 0)
 387       {
 388         vismapping.addShift(npos, contigs[vc]);
 389       }
 390       else
 391       {
 392         // hidden region
 393         vismapping.addShift(npos, contigs[vc] - contigs[vc - 1] + 1);
 394       }
 395       scontigs[vc] = contigs[vc];
 396       scontigs[vc + 1] = contigs[vc + 1];
 397     }
 398
 399     // allocate a roughly sized buffer for the protein sequence
 400     StringBuilder protein = new StringBuilder(seqstring.length() / 2);
 401     String seq = seqstring.replace('U', 'T').replace('u', 'T');
 402     char codon[] = new char[3];
 403     int cdp[] = new int[3], rf = 0, lastnpos = 0, nend;
 404     int aspos = 0;
 405     int resSize = 0;
 406     for (npos = 0, nend = seq.length(); npos < nend; npos++)
 407     {
 408       if (!Comparison.isGap(seq.charAt(npos)))
 409       {
 410         cdp[rf] = npos; // store position
 411         codon[rf++] = seq.charAt(npos); // store base
 412       }
 413       if (rf == 3)
 414       {
 415         /*
 416          * Filled up a reading frame...
 417          */
 418         String aa = ResidueProperties.codonTranslate(new String(codon));
 419         rf = 0;
 420         final String gapString = String.valueOf(gapChar);
 421         if (aa == null)
 422         {
 423           aa = gapString;
 424           if (skipint == null)
 425           {
 426             skipint = new int[]
 427             { cdp[0], cdp[2] };
 428           }
 429           skipint[1] = cdp[2];
 430         }
 431         else
 432         {
 433           if (skipint != null)
 434           {
 435             // edit scontigs
 436             skipint[0] = vismapping.shift(skipint[0]);
 437             skipint[1] = vismapping.shift(skipint[1]);
 438             for (vc = 0; vc < scontigs.length;)
 439             {
 440               if (scontigs[vc + 1] < skipint[0])
 441               {
 442                 // before skipint starts
 443                 vc += 2;
 444                 continue;
 445               }
 446               if (scontigs[vc] > skipint[1])
 447               {
 448                 // finished editing so
 449                 break;
 450               }
 451               // Edit the contig list to include the skipped region which did
 452               // not translate
 453               int[] t;
 454               // from : s1 e1 s2 e2 s3 e3
 455               // to s: s1 e1 s2 k0 k1 e2 s3 e3
 456               // list increases by one unless one boundary (s2==k0 or e2==k1)
 457               // matches, and decreases by one if skipint intersects whole
 458               // visible contig
 459               if (scontigs[vc] <= skipint[0])
 460               {
 461                 if (skipint[0] == scontigs[vc])
 462                 {
 463                   // skipint at start of contig
 464                   // shift the start of this contig
 465                   if (scontigs[vc + 1] > skipint[1])
 466                   {
 467                     scontigs[vc] = skipint[1];
 468                     vc += 2;
 469                   }
 470                   else
 471                   {
 472                     if (scontigs[vc + 1] == skipint[1])
 473                     {
 474                       // remove the contig
 475                       t = new int[scontigs.length - 2];
 476                       if (vc > 0)
 477                       {
 478                         System.arraycopy(scontigs, 0, t, 0, vc - 1);
 479                       }
 480                       if (vc + 2 < t.length)
 481                       {
 482                         System.arraycopy(scontigs, vc + 2, t, vc, t.length
 483                                 - vc + 2);
 484                       }
 485                       scontigs = t;
 486                     }
 487                     else
 488                     {
 489                       // truncate contig to before the skipint region
 490                       scontigs[vc + 1] = skipint[0] - 1;
 491                       vc += 2;
 492                     }
 493                   }
 494                 }
 495                 else
 496                 {
 497                   // scontig starts before start of skipint
 498                   if (scontigs[vc + 1] < skipint[1])
 499                   {
 500                     // skipint truncates end of scontig
 501                     scontigs[vc + 1] = skipint[0] - 1;
 502                     vc += 2;
 503                   }
 504                   else
 505                   {
 506                     // divide region to new contigs
 507                     t = new int[scontigs.length + 2];
 508                     System.arraycopy(scontigs, 0, t, 0, vc + 1);
 509                     t[vc + 1] = skipint[0];
 510                     t[vc + 2] = skipint[1];
 511                     System.arraycopy(scontigs, vc + 1, t, vc + 3,
 512                             scontigs.length - (vc + 1));
 513                     scontigs = t;
 514                     vc += 4;
 515                   }
 516                 }
 517               }
 518             }
 519             skip.add(skipint);
 520             skipint = null;
 521           }
 522           if (aa.equals("STOP"))
 523           {
 524             aa = STOP_X;
 525           }
 526           resSize++;
 527         }
 528         // insert gaps prior to this codon - if necessary
 529         boolean findpos = true;
 530         while (findpos)
 531         {
 532           // expand the codons array if necessary
 533           alignedCodons = checkCodonFrameWidth(alignedCodons, aspos);
 534
 535           /*
 536            * Compare this codon's base positions with those currently aligned to
 537            * this column in the translation.
 538            */
 539           final int compareCodonPos = Dna.compareCodonPos(cdp,
 540                   alignedCodons[aspos]);
 541           // debug
 542           System.out.println(seq + "/" + aa + " codons: "
 543                   + Arrays.deepToString(alignedCodons));
 544           for (SequenceI s : proteinSeqs)
 545           {
 546             System.out.println(s.getSequenceAsString());
 547           }
 548           System.out
 549                   .println(("Compare " + Arrays.toString(cdp) + " at pos "
 550                           + aspos + " with "
 551                           + Arrays.toString(alignedCodons[aspos]) + " got " + compareCodonPos));
 552           // end debug
 553           switch (compareCodonPos)
 554           {
 555           case -1:
 556
 557             /*
 558              * This codon should precede the mapped positions - need to insert a
 559              * gap in all prior sequences.
 560              */
 561             alignedCodons = insertAAGap(aspos, alignedCodons, proteinSeqs);
 562             findpos = false;
 563             break;
 564
 565           case +1:
 566
 567             /*
 568              * This codon belongs after the aligned codons at aspos. Prefix it
 569              * with a gap and try the next position.
 570              */
 571             aa = gapString + aa;
 572             aspos++;
 573             break;
 574
 575           case 0:
 576
 577             /*
 578              * Exact match - codon 'belongs' at this translated position.
 579              */
 580             findpos = false;
 581           }
 582         }
 583         protein.append(aa);
 584         lastnpos = npos;
 585         if (alignedCodons[aspos] == null)
 586         {
 587           // mark this column as aligning to this aligned reading frame
 588           alignedCodons[aspos] = new int[]
 589           { cdp[0], cdp[1], cdp[2] };
 590         }
 591         else if (!Arrays.equals(alignedCodons[aspos], cdp))
 592         {
 593           throw new IllegalStateException("Tried to coalign "
 594                   + Arrays.asList(alignedCodons[aspos], cdp));
 595         }
 596         System.out.println(aspos + "/" + aaWidth);
 597         if (aspos >= aaWidth)
 598         {
 599           // update maximum alignment width
 600           // (we can do this without calling checkCodonFrameWidth because it was
 601           // already done above)
 602           aaWidth = aspos;
 603         }
 604         // ready for next translated reading frame alignment position (if any)
 605         aspos++;
 606       }
 607     }
 608     if (resSize > 0)
 609     {
 610       SequenceI newseq = new Sequence(selection.getName(),
 611               protein.toString());
 612       if (rf != 0)
 613       {
 614         final String errMsg = "trimming contigs for incomplete terminal codon.";
 615         if (Cache.log != null)
 616         {
 617           Cache.log.debug(errMsg);
 618         }
 619         else
 620         {
 621           System.err.println(errMsg);
 622         }
 623         // map and trim contigs to ORF region
 624         vc = scontigs.length - 1;
 625         lastnpos = vismapping.shift(lastnpos); // place npos in context of
 626         // whole dna alignment (rather
 627         // than visible contigs)
 628         // incomplete ORF could be broken over one or two visible contig
 629         // intervals.
 630         while (vc >= 0 && scontigs[vc] > lastnpos)
 631         {
 632           if (vc > 0 && scontigs[vc - 1] > lastnpos)
 633           {
 634             vc -= 2;
 635           }
 636           else
 637           {
 638             // correct last interval in list.
 639             scontigs[vc] = lastnpos;
 640           }
 641         }
 642
 643         if (vc > 0 && (vc + 1) < scontigs.length)
 644         {
 645           // truncate map list to just vc elements
 646           int t[] = new int[vc + 1];
 647           System.arraycopy(scontigs, 0, t, 0, vc + 1);
 648           scontigs = t;
 649         }
 650         if (vc <= 0)
 651         {
 652           scontigs = null;
 653         }
 654       }
 655       if (scontigs != null)
 656       {
 657         npos = 0;
 658         // map scontigs to actual sequence positions on selection
 659         for (vc = 0; vc < scontigs.length; vc += 2)
 660         {
 661           scontigs[vc] = selection.findPosition(scontigs[vc]); // not from 1!
 662           scontigs[vc + 1] = selection.findPosition(scontigs[vc + 1]); // exclusive
 663           if (scontigs[vc + 1] == selection.getEnd())
 664           {
 665             break;
 666           }
 667         }
 668         // trim trailing empty intervals.
 669         if ((vc + 2) < scontigs.length)
 670         {
 671           int t[] = new int[vc + 2];
 672           System.arraycopy(scontigs, 0, t, 0, vc + 2);
 673           scontigs = t;
 674         }
 675         /*
 676          * delete intervals in scontigs which are not translated. 1. map skip
 677          * into sequence position intervals 2. truncate existing ranges and add
 678          * new ranges to exclude untranslated regions. if (skip.size()>0) {
 679          * Vector narange = new Vector(); for (vc=0; vc<scontigs.length; vc++) {
 680          * narange.addElement(new int[] {scontigs[vc]}); } int sint=0,iv[]; vc =
 681          * 0; while (sint<skip.size()) { skipint = (int[]) skip.elementAt(sint);
 682          * do { iv = (int[]) narange.elementAt(vc); if (iv[0]>=skipint[0] &&
 683          * iv[0]<=skipint[1]) { if (iv[0]==skipint[0]) { // delete beginning of
 684          * range } else { // truncate range and create new one if necessary iv =
 685          * (int[]) narange.elementAt(vc+1); if (iv[0]<=skipint[1]) { // truncate
 686          * range iv[0] = skipint[1]; } else { } } } else if (iv[0]<skipint[0]) {
 687          * iv = (int[]) narange.elementAt(vc+1); } } while (iv[0]) } }
 688          */
 689         MapList map = new MapList(scontigs, new int[]
 690         { 1, resSize }, 3, 1);
 691
 692         transferCodedFeatures(selection, newseq, map, null, null);
 693
 694         /*
 695          * Construct a dataset sequence for our new peptide.
 696          */
 697         SequenceI rseq = newseq.deriveSequence();
 698
 699         /*
 700          * Store a mapping (between the dataset sequences for the two
 701          * sequences).
 702          */
 703         // SIDE-EFFECT: acf stores the aligned sequence reseq; to remove!
 704         acf.addMap(selection, rseq, map);
 705         return rseq;
 706       }
 707     }
 708     // register the mapping somehow
 709     //
 710     return null;
 711   }
 712
 713   /**
 714    * Insert a gap into the aligned proteins and the codon mapping array.
 715    *
 716    * @param pos
 717    * @param alignedCodons
 718    * @param proteinSeqs
 719    * @return
 720    */
 721   protected int[][] insertAAGap(int pos, int[][] alignedCodons,
 722           List<SequenceI> proteinSeqs)
 723   {
 724     System.out.println("insertAAGap " + pos + "/" + proteinSeqs.size());
 725     aaWidth++;
 726     for (SequenceI seq : proteinSeqs)
 727     {
 728       seq.insertCharAt(pos, gapChar);
 729     }
 730
 731     int[][] resized = checkCodonFrameWidth(alignedCodons, pos);
 732     if (pos < aaWidth)
 733     {
 734       aaWidth++;
 735       System.arraycopy(resized, pos, resized, pos + 1, resized.length - pos
 736               - 1);
 737       resized[pos] = null; // clear so new codon position can be marked.
 738     }
 739     return resized;
 740   }
 741
 742   /**
 743    * Check the codons array is big enough to accommodate the given position, if
 744    * not resize it.
 745    *
 746    * @param alignedCodons
 747    * @param aspos
 748    * @return the resized array (or the original if no resize needed)
 749    */
 750   protected static int[][] checkCodonFrameWidth(int[][] alignedCodons,
 751           int aspos)
 752   {
 753     // TODO why not codons.length < aspos ?
 754     // should codons expand if length is 2 or 3 and aslen==2 ?
 755     System.out.println("Checking " + alignedCodons.length + "/" + aspos);
 756     if (alignedCodons.length <= aspos + 1)
 757     {
 758       // probably never have to do this ?
 759       int[][] c = new int[alignedCodons.length + 10][];
 760       for (int i = 0; i < alignedCodons.length; i++)
 761       {
 762         c[i] = alignedCodons[i];
 763       }
 764       return c;
 765     }
 766     return alignedCodons;
 767   }
 768
 769   /**
 770    * Given a peptide newly translated from a dna sequence, copy over and set any
 771    * features on the peptide from the DNA. If featureTypes is null, all features
 772    * on the dna sequence are searched (rather than just the displayed ones), and
 773    * similarly for featureGroups.
 774    *
 775    * @param dna
 776    * @param pep
 777    * @param map
 778    * @param featureTypes
 779    *          hash who's keys are the displayed feature type strings
 780    * @param featureGroups
 781    *          hash where keys are feature groups and values are Boolean objects
 782    *          indicating if they are displayed.
 783    */
 784   private static void transferCodedFeatures(SequenceI dna, SequenceI pep,
 785           MapList map, Hashtable featureTypes, Hashtable featureGroups)
 786   {
 787     SequenceFeature[] sf = (dna.getDatasetSequence() != null ? dna
 788             .getDatasetSequence() : dna).getSequenceFeatures();
 789     Boolean fgstate;
 790     DBRefEntry[] dnarefs = DBRefUtils.selectRefs(dna.getDBRef(),
 791             DBRefSource.DNACODINGDBS);
 792     if (dnarefs != null)
 793     {
 794       // intersect with pep
 795       for (int d = 0; d < dnarefs.length; d++)
 796       {
 797         Mapping mp = dnarefs[d].getMap();
 798         if (mp != null)
 799         {
 800         }
 801       }
 802     }
 803     if (sf != null)
 804     {
 805       for (int f = 0; f < sf.length; f++)
 806       {
 807         fgstate = (featureGroups == null) ? null : ((Boolean) featureGroups
 808                 .get(sf[f].featureGroup));
 809         if ((featureTypes == null || featureTypes.containsKey(sf[f]
 810                 .getType())) && (fgstate == null || fgstate.booleanValue()))
 811         {
 812           if (FeatureProperties.isCodingFeature(null, sf[f].getType()))
 813           {
 814             // if (map.intersectsFrom(sf[f].begin, sf[f].end))
 815             {
 816
 817             }
 818           }
 819         }
 820       }
 821     }
 822   }
 823 }