src/jalview/io/FeaturesFile.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.io;
  22
  23 import jalview.analysis.AlignmentUtils;
  24 import jalview.analysis.SequenceIdMatcher;
  25 import jalview.api.AlignViewportI;
  26 import jalview.api.FeatureColourI;
  27 import jalview.api.FeaturesSourceI;
  28 import jalview.datamodel.AlignedCodonFrame;
  29 import jalview.datamodel.Alignment;
  30 import jalview.datamodel.AlignmentI;
  31 import jalview.datamodel.SequenceDummy;
  32 import jalview.datamodel.SequenceFeature;
  33 import jalview.datamodel.SequenceI;
  34 import jalview.io.gff.GffHelperBase;
  35 import jalview.io.gff.GffHelperFactory;
  36 import jalview.io.gff.GffHelperI;
  37 import jalview.schemes.FeatureColour;
  38 import jalview.schemes.UserColourScheme;
  39 import jalview.util.MapList;
  40 import jalview.util.ParseHtmlBodyAndLinks;
  41 import jalview.util.StringUtils;
  42
  43 import java.io.IOException;
  44 import java.util.ArrayList;
  45 import java.util.Arrays;
  46 import java.util.HashMap;
  47 import java.util.Iterator;
  48 import java.util.List;
  49 import java.util.Map;
  50 import java.util.Map.Entry;
  51
  52 /**
  53  * Parses and writes features files, which may be in Jalview, GFF2 or GFF3
  54  * format. These are tab-delimited formats but with differences in the use of
  55  * columns.
  56  *
  57  * A Jalview feature file may define feature colours and then declare that the
  58  * remainder of the file is in GFF format with the line 'GFF'.
  59  *
  60  * GFF3 files may include alignment mappings for features, which Jalview will
  61  * attempt to model, and may include sequence data following a ##FASTA line.
  62  *
  63  *
  64  * @author AMW
  65  * @author jbprocter
  66  * @author gmcarstairs
  67  */
  68 public class FeaturesFile extends AlignFile implements FeaturesSourceI
  69 {
  70   private static final String ID_NOT_SPECIFIED = "ID_NOT_SPECIFIED";
  71
  72   private static final String NOTE = "Note";
  73
  74   protected static final String TAB = "\t";
  75
  76   protected static final String GFF_VERSION = "##gff-version";
  77
  78   private AlignmentI lastmatchedAl = null;
  79
  80   private SequenceIdMatcher matcher = null;
  81
  82   protected AlignmentI dataset;
  83
  84   protected int gffVersion;
  85
  86   /**
  87    * Creates a new FeaturesFile object.
  88    */
  89   public FeaturesFile()
  90   {
  91   }
  92
  93   /**
  94    * Constructor which does not parse the file immediately
  95    *
  96    * @param inFile
  97    * @param paste
  98    * @throws IOException
  99    */
 100   public FeaturesFile(String inFile, DataSourceType paste)
 101           throws IOException
 102   {
 103     super(false, inFile, paste);
 104   }
 105
 106   /**
 107    * @param source
 108    * @throws IOException
 109    */
 110   public FeaturesFile(FileParse source) throws IOException
 111   {
 112     super(source);
 113   }
 114
 115   /**
 116    * Constructor that optionally parses the file immediately
 117    *
 118    * @param parseImmediately
 119    * @param inFile
 120    * @param type
 121    * @throws IOException
 122    */
 123   public FeaturesFile(boolean parseImmediately, String inFile,
 124           DataSourceType type)
 125           throws IOException
 126   {
 127     super(parseImmediately, inFile, type);
 128   }
 129
 130   /**
 131    * Parse GFF or sequence features file using case-independent matching,
 132    * discarding URLs
 133    *
 134    * @param align
 135    *          - alignment/dataset containing sequences that are to be annotated
 136    * @param colours
 137    *          - hashtable to store feature colour definitions
 138    * @param removeHTML
 139    *          - process html strings into plain text
 140    * @return true if features were added
 141    */
 142   public boolean parse(AlignmentI align,
 143           Map<String, FeatureColourI> colours, boolean removeHTML)
 144   {
 145     return parse(align, colours, removeHTML, false);
 146   }
 147
 148   /**
 149    * Extends the default addProperties by also adding peptide-to-cDNA mappings
 150    * (if any) derived while parsing a GFF file
 151    */
 152   @Override
 153   public void addProperties(AlignmentI al)
 154   {
 155     super.addProperties(al);
 156     if (dataset != null && dataset.getCodonFrames() != null)
 157     {
 158       AlignmentI ds = (al.getDataset() == null) ? al : al.getDataset();
 159       for (AlignedCodonFrame codons : dataset.getCodonFrames())
 160       {
 161         ds.addCodonFrame(codons);
 162       }
 163     }
 164   }
 165
 166   /**
 167    * Parse GFF or Jalview format sequence features file
 168    *
 169    * @param align
 170    *          - alignment/dataset containing sequences that are to be annotated
 171    * @param colours
 172    *          - hashtable to store feature colour definitions
 173    * @param removeHTML
 174    *          - process html strings into plain text
 175    * @param relaxedIdmatching
 176    *          - when true, ID matches to compound sequence IDs are allowed
 177    * @return true if features were added
 178    */
 179   public boolean parse(AlignmentI align,
 180           Map<String, FeatureColourI> colours, boolean removeHTML,
 181           boolean relaxedIdmatching)
 182   {
 183     Map<String, String> gffProps = new HashMap<String, String>();
 184     /*
 185      * keep track of any sequences we try to create from the data
 186      */
 187     List<SequenceI> newseqs = new ArrayList<SequenceI>();
 188
 189     String line = null;
 190     try
 191     {
 192       String[] gffColumns;
 193       String featureGroup = null;
 194
 195       while ((line = nextLine()) != null)
 196       {
 197         // skip comments/process pragmas
 198         if (line.length() == 0 || line.startsWith("#"))
 199         {
 200           if (line.toLowerCase().startsWith("##"))
 201           {
 202             processGffPragma(line, gffProps, align, newseqs);
 203           }
 204           continue;
 205         }
 206
 207         gffColumns = line.split("\\t"); // tab as regex
 208         if (gffColumns.length == 1)
 209         {
 210           if (line.trim().equalsIgnoreCase("GFF"))
 211           {
 212             /*
 213              * Jalview features file with appended GFF
 214              * assume GFF2 (though it may declare ##gff-version 3)
 215              */
 216             gffVersion = 2;
 217             continue;
 218           }
 219         }
 220
 221         if (gffColumns.length > 1 && gffColumns.length < 4)
 222         {
 223           /*
 224            * if 2 or 3 tokens, we anticipate either 'startgroup', 'endgroup' or
 225            * a feature type colour specification
 226            */
 227           String ft = gffColumns[0];
 228           if (ft.equalsIgnoreCase("startgroup"))
 229           {
 230             featureGroup = gffColumns[1];
 231           }
 232           else if (ft.equalsIgnoreCase("endgroup"))
 233           {
 234             // We should check whether this is the current group,
 235             // but at present there's no way of showing more than 1 group
 236             featureGroup = null;
 237           }
 238           else
 239           {
 240             String colscheme = gffColumns[1];
 241             FeatureColourI colour = FeatureColour
 242                     .parseJalviewFeatureColour(colscheme);
 243             if (colour != null)
 244             {
 245               colours.put(ft, colour);
 246             }
 247           }
 248           continue;
 249         }
 250
 251         /*
 252          * if not a comment, GFF pragma, startgroup, endgroup or feature
 253          * colour specification, that just leaves a feature details line
 254          * in either Jalview or GFF format
 255          */
 256         if (gffVersion == 0)
 257         {
 258           parseJalviewFeature(line, gffColumns, align, colours, removeHTML,
 259                   relaxedIdmatching, featureGroup);
 260         }
 261         else
 262         {
 263           parseGff(gffColumns, align, relaxedIdmatching, newseqs);
 264         }
 265       }
 266       resetMatcher();
 267     } catch (Exception ex)
 268     {
 269       // should report somewhere useful for UI if necessary
 270       warningMessage = ((warningMessage == null) ? "" : warningMessage)
 271               + "Parsing error at\n" + line;
 272       System.out.println("Error parsing feature file: " + ex + "\n" + line);
 273       ex.printStackTrace(System.err);
 274       resetMatcher();
 275       return false;
 276     }
 277
 278     /*
 279      * experimental - add any dummy sequences with features to the alignment
 280      * - we need them for Ensembl feature extraction - though maybe not otherwise
 281      */
 282     for (SequenceI newseq : newseqs)
 283     {
 284       if (newseq.getSequenceFeatures() != null)
 285       {
 286         align.addSequence(newseq);
 287       }
 288     }
 289     return true;
 290   }
 291
 292   /**
 293    * Try to parse a Jalview format feature specification and add it as a
 294    * sequence feature to any matching sequences in the alignment. Returns true
 295    * if successful (a feature was added), or false if not.
 296    *
 297    * @param line
 298    * @param gffColumns
 299    * @param alignment
 300    * @param featureColours
 301    * @param removeHTML
 302    * @param relaxedIdmatching
 303    * @param featureGroup
 304    */
 305   protected boolean parseJalviewFeature(String line, String[] gffColumns,
 306           AlignmentI alignment, Map<String, FeatureColourI> featureColours,
 307           boolean removeHTML, boolean relaxedIdMatching, String featureGroup)
 308   {
 309     /*
 310      * tokens: description seqid seqIndex start end type [score]
 311      */
 312     if (gffColumns.length < 6)
 313     {
 314       System.err.println("Ignoring feature line '" + line
 315               + "' with too few columns (" + gffColumns.length + ")");
 316       return false;
 317     }
 318     String desc = gffColumns[0];
 319     String seqId = gffColumns[1];
 320     SequenceI seq = findSequence(seqId, alignment, null, relaxedIdMatching);
 321
 322     if (!ID_NOT_SPECIFIED.equals(seqId))
 323     {
 324       seq = findSequence(seqId, alignment, null, relaxedIdMatching);
 325     }
 326     else
 327     {
 328       seqId = null;
 329       seq = null;
 330       String seqIndex = gffColumns[2];
 331       try
 332       {
 333         int idx = Integer.parseInt(seqIndex);
 334         seq = alignment.getSequenceAt(idx);
 335       } catch (NumberFormatException ex)
 336       {
 337         System.err.println("Invalid sequence index: " + seqIndex);
 338       }
 339     }
 340
 341     if (seq == null)
 342     {
 343       System.out.println("Sequence not found: " + line);
 344       return false;
 345     }
 346
 347     int startPos = Integer.parseInt(gffColumns[3]);
 348     int endPos = Integer.parseInt(gffColumns[4]);
 349
 350     String ft = gffColumns[5];
 351
 352     if (!featureColours.containsKey(ft))
 353     {
 354       /*
 355        * Perhaps an old style groups file with no colours -
 356        * synthesize a colour from the feature type
 357        */
 358       UserColourScheme ucs = new UserColourScheme(ft);
 359       featureColours.put(ft, new FeatureColour(ucs.findColour('A')));
 360     }
 361     SequenceFeature sf = new SequenceFeature(ft, desc, "", startPos,
 362             endPos, featureGroup);
 363     if (gffColumns.length > 6)
 364     {
 365       float score = Float.NaN;
 366       try
 367       {
 368         score = new Float(gffColumns[6]).floatValue();
 369         // update colourgradient bounds if allowed to
 370       } catch (NumberFormatException ex)
 371       {
 372         // leave as NaN
 373       }
 374       sf.setScore(score);
 375     }
 376
 377     parseDescriptionHTML(sf, removeHTML);
 378
 379     seq.addSequenceFeature(sf);
 380
 381     while (seqId != null
 382             && (seq = alignment.findName(seq, seqId, false)) != null)
 383     {
 384       seq.addSequenceFeature(new SequenceFeature(sf));
 385     }
 386     return true;
 387   }
 388
 389   /**
 390    * clear any temporary handles used to speed up ID matching
 391    */
 392   protected void resetMatcher()
 393   {
 394     lastmatchedAl = null;
 395     matcher = null;
 396   }
 397
 398   /**
 399    * Returns a sequence matching the given id, as follows
 400    * <ul>
 401    * <li>strict matching is on exact sequence name</li>
 402    * <li>relaxed matching allows matching on a token within the sequence name,
 403    * or a dbxref</li>
 404    * <li>first tries to find a match in the alignment sequences</li>
 405    * <li>else tries to find a match in the new sequences already generated while
 406    * parsing the features file</li>
 407    * <li>else creates a new placeholder sequence, adds it to the new sequences
 408    * list, and returns it</li>
 409    * </ul>
 410    *
 411    * @param seqId
 412    * @param align
 413    * @param newseqs
 414    * @param relaxedIdMatching
 415    *
 416    * @return
 417    */
 418   protected SequenceI findSequence(String seqId, AlignmentI align,
 419           List<SequenceI> newseqs, boolean relaxedIdMatching)
 420   {
 421     // TODO encapsulate in SequenceIdMatcher, share the matcher
 422     // with the GffHelper (removing code duplication)
 423     SequenceI match = null;
 424     if (relaxedIdMatching)
 425     {
 426       if (lastmatchedAl != align)
 427       {
 428         lastmatchedAl = align;
 429         matcher = new SequenceIdMatcher(align.getSequencesArray());
 430         if (newseqs != null)
 431         {
 432           matcher.addAll(newseqs);
 433         }
 434       }
 435       match = matcher.findIdMatch(seqId);
 436     }
 437     else
 438     {
 439       match = align.findName(seqId, true);
 440       if (match == null && newseqs != null)
 441       {
 442         for (SequenceI m : newseqs)
 443         {
 444           if (seqId.equals(m.getName()))
 445           {
 446             return m;
 447           }
 448         }
 449       }
 450
 451     }
 452     if (match == null && newseqs != null)
 453     {
 454       match = new SequenceDummy(seqId);
 455       if (relaxedIdMatching)
 456       {
 457         matcher.addAll(Arrays.asList(new SequenceI[] { match }));
 458       }
 459       // add dummy sequence to the newseqs list
 460       newseqs.add(match);
 461     }
 462     return match;
 463   }
 464
 465   public void parseDescriptionHTML(SequenceFeature sf, boolean removeHTML)
 466   {
 467     if (sf.getDescription() == null)
 468     {
 469       return;
 470     }
 471     ParseHtmlBodyAndLinks parsed = new ParseHtmlBodyAndLinks(
 472             sf.getDescription(), removeHTML, newline);
 473
 474     sf.description = (removeHTML) ? parsed.getNonHtmlContent()
 475             : sf.description;
 476     for (String link : parsed.getLinks())
 477     {
 478       sf.addLink(link);
 479     }
 480
 481   }
 482
 483   /**
 484    * generate a features file for seqs includes non-pos features by default.
 485    *
 486    * @param sequences
 487    *          source of sequence features
 488    * @param visible
 489    *          hash of feature types and colours
 490    * @return features file contents
 491    */
 492   public String printJalviewFormat(SequenceI[] sequences,
 493           Map<String, FeatureColourI> visible)
 494   {
 495     return printJalviewFormat(sequences, visible, true, true);
 496   }
 497
 498   /**
 499    * generate a features file for seqs with colours from visible (if any)
 500    *
 501    * @param sequences
 502    *          source of features
 503    * @param visible
 504    *          hash of Colours for each feature type
 505    * @param visOnly
 506    *          when true only feature types in 'visible' will be output
 507    * @param nonpos
 508    *          indicates if non-positional features should be output (regardless
 509    *          of group or type)
 510    * @return features file contents
 511    */
 512   public String printJalviewFormat(SequenceI[] sequences,
 513           Map<String, FeatureColourI> visible, boolean visOnly,
 514           boolean nonpos)
 515   {
 516     StringBuilder out = new StringBuilder(256);
 517     boolean featuresGen = false;
 518     if (visOnly && !nonpos && (visible == null || visible.size() < 1))
 519     {
 520       // no point continuing.
 521       return "No Features Visible";
 522     }
 523
 524     if (visible != null && visOnly)
 525     {
 526       // write feature colours only if we're given them and we are generating
 527       // viewed features
 528       // TODO: decide if feature links should also be written here ?
 529       Iterator<String> en = visible.keySet().iterator();
 530       while (en.hasNext())
 531       {
 532         String featureType = en.next().toString();
 533         FeatureColourI colour = visible.get(featureType);
 534         out.append(colour.toJalviewFormat(featureType)).append(newline);
 535       }
 536     }
 537
 538     // Work out which groups are both present and visible
 539     List<String> groups = new ArrayList<String>();
 540     int groupIndex = 0;
 541     boolean isnonpos = false;
 542
 543     SequenceFeature[] features;
 544     for (int i = 0; i < sequences.length; i++)
 545     {
 546       features = sequences[i].getSequenceFeatures();
 547       if (features != null)
 548       {
 549         for (int j = 0; j < features.length; j++)
 550         {
 551           isnonpos = features[j].begin == 0 && features[j].end == 0;
 552           if ((!nonpos && isnonpos)
 553                   || (!isnonpos && visOnly && !visible
 554                           .containsKey(features[j].type)))
 555           {
 556             continue;
 557           }
 558
 559           if (features[j].featureGroup != null
 560                   && !groups.contains(features[j].featureGroup))
 561           {
 562             groups.add(features[j].featureGroup);
 563           }
 564         }
 565       }
 566     }
 567
 568     String group = null;
 569     do
 570     {
 571       if (groups.size() > 0 && groupIndex < groups.size())
 572       {
 573         group = groups.get(groupIndex);
 574         out.append(newline);
 575         out.append("STARTGROUP").append(TAB);
 576         out.append(group);
 577         out.append(newline);
 578       }
 579       else
 580       {
 581         group = null;
 582       }
 583
 584       for (int i = 0; i < sequences.length; i++)
 585       {
 586         features = sequences[i].getSequenceFeatures();
 587         if (features != null)
 588         {
 589           for (SequenceFeature sequenceFeature : features)
 590           {
 591             isnonpos = sequenceFeature.begin == 0
 592                     && sequenceFeature.end == 0;
 593             if ((!nonpos && isnonpos)
 594                     || (!isnonpos && visOnly && !visible
 595                             .containsKey(sequenceFeature.type)))
 596             {
 597               // skip if feature is nonpos and we ignore them or if we only
 598               // output visible and it isn't non-pos and it's not visible
 599               continue;
 600             }
 601
 602             if (group != null
 603                     && (sequenceFeature.featureGroup == null || !sequenceFeature.featureGroup
 604                             .equals(group)))
 605             {
 606               continue;
 607             }
 608
 609             if (group == null && sequenceFeature.featureGroup != null)
 610             {
 611               continue;
 612             }
 613             // we have features to output
 614             featuresGen = true;
 615             if (sequenceFeature.description == null
 616                     || sequenceFeature.description.equals(""))
 617             {
 618               out.append(sequenceFeature.type).append(TAB);
 619             }
 620             else
 621             {
 622               if (sequenceFeature.links != null
 623                       && sequenceFeature.getDescription().indexOf("<html>") == -1)
 624               {
 625                 out.append("<html>");
 626               }
 627
 628               out.append(sequenceFeature.description);
 629               if (sequenceFeature.links != null)
 630               {
 631                 for (int l = 0; l < sequenceFeature.links.size(); l++)
 632                 {
 633                   String label = sequenceFeature.links.elementAt(l);
 634                   String href = label.substring(label.indexOf("|") + 1);
 635                   label = label.substring(0, label.indexOf("|"));
 636
 637                   if (sequenceFeature.description.indexOf(href) == -1)
 638                   {
 639                     out.append(" <a href=\"" + href + "\">" + label
 640                             + "</a>");
 641                   }
 642                 }
 643
 644                 if (sequenceFeature.getDescription().indexOf("</html>") == -1)
 645                 {
 646                   out.append("</html>");
 647                 }
 648               }
 649
 650               out.append(TAB);
 651             }
 652             out.append(sequences[i].getName());
 653             out.append("\t-1\t");
 654             out.append(sequenceFeature.begin);
 655             out.append(TAB);
 656             out.append(sequenceFeature.end);
 657             out.append(TAB);
 658             out.append(sequenceFeature.type);
 659             if (!Float.isNaN(sequenceFeature.score))
 660             {
 661               out.append(TAB);
 662               out.append(sequenceFeature.score);
 663             }
 664             out.append(newline);
 665           }
 666         }
 667       }
 668
 669       if (group != null)
 670       {
 671         out.append("ENDGROUP").append(TAB);
 672         out.append(group);
 673         out.append(newline);
 674         groupIndex++;
 675       }
 676       else
 677       {
 678         break;
 679       }
 680
 681     } while (groupIndex < groups.size() + 1);
 682
 683     if (!featuresGen)
 684     {
 685       return "No Features Visible";
 686     }
 687
 688     return out.toString();
 689   }
 690
 691   /**
 692    * Parse method that is called when a GFF file is dragged to the desktop
 693    */
 694   @Override
 695   public void parse()
 696   {
 697     AlignViewportI av = getViewport();
 698     if (av != null)
 699     {
 700       if (av.getAlignment() != null)
 701       {
 702         dataset = av.getAlignment().getDataset();
 703       }
 704       if (dataset == null)
 705       {
 706         // working in the applet context ?
 707         dataset = av.getAlignment();
 708       }
 709     }
 710     else
 711     {
 712       dataset = new Alignment(new SequenceI[] {});
 713     }
 714
 715     Map<String, FeatureColourI> featureColours = new HashMap<String, FeatureColourI>();
 716     boolean parseResult = parse(dataset, featureColours, false, true);
 717     if (!parseResult)
 718     {
 719       // pass error up somehow
 720     }
 721     if (av != null)
 722     {
 723       // update viewport with the dataset data ?
 724     }
 725     else
 726     {
 727       setSeqs(dataset.getSequencesArray());
 728     }
 729   }
 730
 731   /**
 732    * Implementation of unused abstract method
 733    *
 734    * @return error message
 735    */
 736   @Override
 737   public String print(SequenceI[] sqs, boolean jvsuffix)
 738   {
 739     System.out.println("Use printGffFormat() or printJalviewFormat()");
 740     return null;
 741   }
 742
 743   /**
 744    * Returns features output in GFF2 format, including hidden and non-positional
 745    * features
 746    *
 747    * @param sequences
 748    *          the sequences whose features are to be output
 749    * @param visible
 750    *          a map whose keys are the type names of visible features
 751    * @return
 752    */
 753   public String printGffFormat(SequenceI[] sequences,
 754           Map<String, FeatureColourI> visible)
 755   {
 756     return printGffFormat(sequences, visible, true, true);
 757   }
 758
 759   /**
 760    * Returns features output in GFF2 format
 761    *
 762    * @param sequences
 763    *          the sequences whose features are to be output
 764    * @param visible
 765    *          a map whose keys are the type names of visible features
 766    * @param outputVisibleOnly
 767    * @param includeNonPositionalFeatures
 768    * @return
 769    */
 770   public String printGffFormat(SequenceI[] sequences,
 771           Map<String, FeatureColourI> visible, boolean outputVisibleOnly,
 772           boolean includeNonPositionalFeatures)
 773   {
 774     StringBuilder out = new StringBuilder(256);
 775     int version = gffVersion == 0 ? 2 : gffVersion;
 776     out.append(String.format("%s %d\n", GFF_VERSION, version));
 777     String source;
 778     boolean isnonpos;
 779     for (SequenceI seq : sequences)
 780     {
 781       SequenceFeature[] features = seq.getSequenceFeatures();
 782       if (features != null)
 783       {
 784         for (SequenceFeature sf : features)
 785         {
 786           isnonpos = sf.begin == 0 && sf.end == 0;
 787           if (!includeNonPositionalFeatures && isnonpos)
 788           {
 789             /*
 790              * ignore non-positional features if not wanted
 791              */
 792             continue;
 793           }
 794           // TODO why the test !isnonpos here?
 795           // what about not visible non-positional features?
 796           if (!isnonpos && outputVisibleOnly
 797                   && !visible.containsKey(sf.type))
 798           {
 799             /*
 800              * ignore not visible features if not wanted
 801              */
 802             continue;
 803           }
 804
 805           source = sf.featureGroup;
 806           if (source == null)
 807           {
 808             source = sf.getDescription();
 809           }
 810
 811           out.append(seq.getName());
 812           out.append(TAB);
 813           out.append(source);
 814           out.append(TAB);
 815           out.append(sf.type);
 816           out.append(TAB);
 817           out.append(sf.begin);
 818           out.append(TAB);
 819           out.append(sf.end);
 820           out.append(TAB);
 821           out.append(sf.score);
 822           out.append(TAB);
 823
 824           int strand = sf.getStrand();
 825           out.append(strand == 1 ? "+" : (strand == -1 ? "-" : "."));
 826           out.append(TAB);
 827
 828           String phase = sf.getPhase();
 829           out.append(phase == null ? "." : phase);
 830
 831           // miscellaneous key-values (GFF column 9)
 832           String attributes = sf.getAttributes();
 833           if (attributes != null)
 834           {
 835             out.append(TAB).append(attributes);
 836           }
 837
 838           out.append(newline);
 839         }
 840       }
 841     }
 842
 843     return out.toString();
 844   }
 845
 846   /**
 847    * Returns a mapping given list of one or more Align descriptors (exonerate
 848    * format)
 849    *
 850    * @param alignedRegions
 851    *          a list of "Align fromStart toStart fromCount"
 852    * @param mapIsFromCdna
 853    *          if true, 'from' is dna, else 'from' is protein
 854    * @param strand
 855    *          either 1 (forward) or -1 (reverse)
 856    * @return
 857    * @throws IOException
 858    */
 859   protected MapList constructCodonMappingFromAlign(
 860           List<String> alignedRegions, boolean mapIsFromCdna, int strand)
 861           throws IOException
 862   {
 863     if (strand == 0)
 864     {
 865       throw new IOException(
 866               "Invalid strand for a codon mapping (cannot be 0)");
 867     }
 868     int regions = alignedRegions.size();
 869     // arrays to hold [start, end] for each aligned region
 870     int[] fromRanges = new int[regions * 2]; // from dna
 871     int[] toRanges = new int[regions * 2]; // to protein
 872     int fromRangesIndex = 0;
 873     int toRangesIndex = 0;
 874
 875     for (String range : alignedRegions)
 876     {
 877       /*
 878        * Align mapFromStart mapToStart mapFromCount
 879        * e.g. if mapIsFromCdna
 880        *     Align 11270 143 120
 881        * means:
 882        *     120 bases from pos 11270 align to pos 143 in peptide
 883        * if !mapIsFromCdna this would instead be
 884        *     Align 143 11270 40
 885        */
 886       String[] tokens = range.split(" ");
 887       if (tokens.length != 3)
 888       {
 889         throw new IOException("Wrong number of fields for Align");
 890       }
 891       int fromStart = 0;
 892       int toStart = 0;
 893       int fromCount = 0;
 894       try
 895       {
 896         fromStart = Integer.parseInt(tokens[0]);
 897         toStart = Integer.parseInt(tokens[1]);
 898         fromCount = Integer.parseInt(tokens[2]);
 899       } catch (NumberFormatException nfe)
 900       {
 901         throw new IOException("Invalid number in Align field: "
 902                 + nfe.getMessage());
 903       }
 904
 905       /*
 906        * Jalview always models from dna to protein, so adjust values if the
 907        * GFF mapping is from protein to dna
 908        */
 909       if (!mapIsFromCdna)
 910       {
 911         fromCount *= 3;
 912         int temp = fromStart;
 913         fromStart = toStart;
 914         toStart = temp;
 915       }
 916       fromRanges[fromRangesIndex++] = fromStart;
 917       fromRanges[fromRangesIndex++] = fromStart + strand * (fromCount - 1);
 918
 919       /*
 920        * If a codon has an intron gap, there will be contiguous 'toRanges';
 921        * this is handled for us by the MapList constructor.
 922        * (It is not clear that exonerate ever generates this case)
 923        */
 924       toRanges[toRangesIndex++] = toStart;
 925       toRanges[toRangesIndex++] = toStart + (fromCount - 1) / 3;
 926     }
 927
 928     return new MapList(fromRanges, toRanges, 3, 1);
 929   }
 930
 931   /**
 932    * Parse a GFF format feature. This may include creating a 'dummy' sequence to
 933    * hold the feature, or for its mapped sequence, or both, to be resolved
 934    * either later in the GFF file (##FASTA section), or when the user loads
 935    * additional sequences.
 936    *
 937    * @param gffColumns
 938    * @param alignment
 939    * @param relaxedIdMatching
 940    * @param newseqs
 941    * @return
 942    */
 943   protected SequenceI parseGff(String[] gffColumns, AlignmentI alignment,
 944           boolean relaxedIdMatching, List<SequenceI> newseqs)
 945   {
 946     /*
 947      * GFF: seqid source type start end score strand phase [attributes]
 948      */
 949     if (gffColumns.length < 5)
 950     {
 951       System.err.println("Ignoring GFF feature line with too few columns ("
 952               + gffColumns.length + ")");
 953       return null;
 954     }
 955
 956     /*
 957      * locate referenced sequence in alignment _or_
 958      * as a forward or external reference (SequenceDummy)
 959      */
 960     String seqId = gffColumns[0];
 961     SequenceI seq = findSequence(seqId, alignment, newseqs,
 962             relaxedIdMatching);
 963
 964     SequenceFeature sf = null;
 965     GffHelperI helper = GffHelperFactory.getHelper(gffColumns);
 966     if (helper != null)
 967     {
 968       try
 969       {
 970         sf = helper.processGff(seq, gffColumns, alignment, newseqs,
 971                 relaxedIdMatching);
 972         if (sf != null)
 973         {
 974           seq.addSequenceFeature(sf);
 975           while ((seq = alignment.findName(seq, seqId, true)) != null)
 976           {
 977             seq.addSequenceFeature(new SequenceFeature(sf));
 978           }
 979         }
 980       } catch (IOException e)
 981       {
 982         System.err.println("GFF parsing failed with: " + e.getMessage());
 983         return null;
 984       }
 985     }
 986
 987     return seq;
 988   }
 989
 990   /**
 991    * Process the 'column 9' data of the GFF file. This is less formally defined,
 992    * and its interpretation will vary depending on the tool that has generated
 993    * it.
 994    *
 995    * @param attributes
 996    * @param sf
 997    */
 998   protected void processGffColumnNine(String attributes, SequenceFeature sf)
 999   {
1000     sf.setAttributes(attributes);
1001
1002     /*
1003      * Parse attributes in column 9 and add them to the sequence feature's
1004      * 'otherData' table; use Note as a best proxy for description
1005      */
1006     char nameValueSeparator = gffVersion == 3 ? '=' : ' ';
1007     // TODO check we don't break GFF2 values which include commas here
1008     Map<String, List<String>> nameValues = GffHelperBase
1009             .parseNameValuePairs(attributes, ";", nameValueSeparator, ",");
1010     for (Entry<String, List<String>> attr : nameValues.entrySet())
1011     {
1012       String values = StringUtils.listToDelimitedString(attr.getValue(),
1013               "; ");
1014       sf.setValue(attr.getKey(), values);
1015       if (NOTE.equals(attr.getKey()))
1016       {
1017         sf.setDescription(values);
1018       }
1019     }
1020   }
1021
1022   /**
1023    * After encountering ##fasta in a GFF3 file, process the remainder of the
1024    * file as FAST sequence data. Any placeholder sequences created during
1025    * feature parsing are updated with the actual sequences.
1026    *
1027    * @param align
1028    * @param newseqs
1029    * @throws IOException
1030    */
1031   protected void processAsFasta(AlignmentI align, List<SequenceI> newseqs)
1032           throws IOException
1033   {
1034     try
1035     {
1036       mark();
1037     } catch (IOException q)
1038     {
1039     }
1040     FastaFile parser = new FastaFile(this);
1041     List<SequenceI> includedseqs = parser.getSeqs();
1042
1043     SequenceIdMatcher smatcher = new SequenceIdMatcher(newseqs);
1044
1045     /*
1046      * iterate over includedseqs, and replacing matching ones with newseqs
1047      * sequences. Generic iterator not used here because we modify
1048      * includedseqs as we go
1049      */
1050     for (int p = 0, pSize = includedseqs.size(); p < pSize; p++)
1051     {
1052       // search for any dummy seqs that this sequence can be used to update
1053       SequenceI includedSeq = includedseqs.get(p);
1054       SequenceI dummyseq = smatcher.findIdMatch(includedSeq);
1055       if (dummyseq != null && dummyseq instanceof SequenceDummy)
1056       {
1057         // probably have the pattern wrong
1058         // idea is that a flyweight proxy for a sequence ID can be created for
1059         // 1. stable reference creation
1060         // 2. addition of annotation
1061         // 3. future replacement by a real sequence
1062         // current pattern is to create SequenceDummy objects - a convenience
1063         // constructor for a Sequence.
1064         // problem is that when promoted to a real sequence, all references
1065         // need to be updated somehow. We avoid that by keeping the same object.
1066         ((SequenceDummy) dummyseq).become(includedSeq);
1067         dummyseq.createDatasetSequence();
1068
1069         /*
1070          * Update mappings so they are now to the dataset sequence
1071          */
1072         for (AlignedCodonFrame mapping : align.getCodonFrames())
1073         {
1074           mapping.updateToDataset(dummyseq);
1075         }
1076
1077         /*
1078          * replace parsed sequence with the realised forward reference
1079          */
1080         includedseqs.set(p, dummyseq);
1081
1082         /*
1083          * and remove from the newseqs list
1084          */
1085         newseqs.remove(dummyseq);
1086       }
1087     }
1088
1089     /*
1090      * finally add sequences to the dataset
1091      */
1092     for (SequenceI seq : includedseqs)
1093     {
1094       // experimental: mapping-based 'alignment' to query sequence
1095       AlignmentUtils.alignSequenceAs(seq, align,
1096               String.valueOf(align.getGapCharacter()), false, true);
1097
1098       // rename sequences if GFF handler requested this
1099       // TODO a more elegant way e.g. gffHelper.postProcess(newseqs) ?
1100       SequenceFeature[] sfs = seq.getSequenceFeatures();
1101       if (sfs != null)
1102       {
1103         String newName = (String) sfs[0].getValue(GffHelperI.RENAME_TOKEN);
1104         if (newName != null)
1105         {
1106           seq.setName(newName);
1107         }
1108       }
1109       align.addSequence(seq);
1110     }
1111   }
1112
1113   /**
1114    * Process a ## directive
1115    *
1116    * @param line
1117    * @param gffProps
1118    * @param align
1119    * @param newseqs
1120    * @throws IOException
1121    */
1122   protected void processGffPragma(String line,
1123           Map<String, String> gffProps, AlignmentI align,
1124           List<SequenceI> newseqs) throws IOException
1125   {
1126     line = line.trim();
1127     if ("###".equals(line))
1128     {
1129       // close off any open 'forward references'
1130       return;
1131     }
1132
1133     String[] tokens = line.substring(2).split(" ");
1134     String pragma = tokens[0];
1135     String value = tokens.length == 1 ? null : tokens[1];
1136
1137     if ("gff-version".equalsIgnoreCase(pragma))
1138     {
1139       if (value != null)
1140       {
1141         try
1142         {
1143           // value may be e.g. "3.1.2"
1144           gffVersion = Integer.parseInt(value.split("\\.")[0]);
1145         } catch (NumberFormatException e)
1146         {
1147           // ignore
1148         }
1149       }
1150     }
1151     else if ("sequence-region".equalsIgnoreCase(pragma))
1152     {
1153       // could capture <seqid start end> if wanted here
1154     }
1155     else if ("feature-ontology".equalsIgnoreCase(pragma))
1156     {
1157       // should resolve against the specified feature ontology URI
1158     }
1159     else if ("attribute-ontology".equalsIgnoreCase(pragma))
1160     {
1161       // URI of attribute ontology - not currently used in GFF3
1162     }
1163     else if ("source-ontology".equalsIgnoreCase(pragma))
1164     {
1165       // URI of source ontology - not currently used in GFF3
1166     }
1167     else if ("species-build".equalsIgnoreCase(pragma))
1168     {
1169       // save URI of specific NCBI taxon version of annotations
1170       gffProps.put("species-build", value);
1171     }
1172     else if ("fasta".equalsIgnoreCase(pragma))
1173     {
1174       // process the rest of the file as a fasta file and replace any dummy
1175       // sequence IDs
1176       processAsFasta(align, newseqs);
1177     }
1178     else
1179     {
1180       System.err.println("Ignoring unknown pragma: " + line);
1181     }
1182   }
1183 }