JAL-2110 field 'cf' moved to stack
[jalview.git] / src / jalview / analysis / CrossRef.java
index d96ab58..288d60e 100644 (file)
@@ -28,7 +28,6 @@ import jalview.datamodel.Mapping;
 import jalview.datamodel.Sequence;
 import jalview.datamodel.SequenceFeature;
 import jalview.datamodel.SequenceI;
-import jalview.util.Comparison;
 import jalview.util.DBRefUtils;
 import jalview.util.MapList;
 import jalview.ws.SequenceFetcherFactory;
@@ -59,6 +58,16 @@ public class CrossRef
   private SequenceI[] fromSeqs;
 
   /**
+   * matcher built from dataset
+   */
+  SequenceIdMatcher matcher;
+
+  /**
+   * sequences found by cross-ref searches to fromSeqs
+   */
+  List<SequenceI> rseqs;
+
+  /**
    * Constructor
    * 
    * @param seqs
@@ -111,6 +120,9 @@ public class CrossRef
    * 
    * @param seq
    *          the sequence whose dbrefs we are searching against
+   * @param fromDna
+   *          when true, context is DNA - so sources identifying protein
+   *          products will be returned.
    * @param sources
    *          a list of sources to add matches to
    */
@@ -128,18 +140,18 @@ public class CrossRef
        * find sequence's direct (dna-to-dna, peptide-to-peptide) xrefs
        */
       DBRefEntry[] lrfs = DBRefUtils.selectDbRefs(fromDna, seq.getDBRefs());
-      List<SequenceI> rseqs = new ArrayList<SequenceI>();
+      List<SequenceI> foundSeqs = new ArrayList<SequenceI>();
 
       /*
        * find sequences in the alignment which xref one of these DBRefs
        * i.e. is xref-ed to a common sequence identifier
        */
-      searchDatasetXrefs(fromDna, seq, lrfs, rseqs, null);
+      searchDatasetXrefs(fromDna, seq, lrfs, foundSeqs, null);
 
       /*
        * add those sequences' (dna-to-peptide or peptide-to-dna) dbref sources
        */
-      for (SequenceI rs : rseqs)
+      for (SequenceI rs : foundSeqs)
       {
         DBRefEntry[] xrs = DBRefUtils
                 .selectDbRefs(!fromDna, rs.getDBRefs());
@@ -195,9 +207,9 @@ public class CrossRef
   public Alignment findXrefSequences(String source, boolean fromDna)
   {
 
-    List<SequenceI> rseqs = new ArrayList<SequenceI>();
+    rseqs = new ArrayList<SequenceI>();
     AlignedCodonFrame cf = new AlignedCodonFrame();
-    SequenceIdMatcher matcher = new SequenceIdMatcher(
+    matcher = new SequenceIdMatcher(
             dataset.getSequences());
 
     for (SequenceI seq : fromSeqs)
@@ -278,12 +290,14 @@ public class CrossRef
               if (fromDna)
               {
                 // map is from dna seq to a protein product
-                cf.addMap(dss, rsq, xref.getMap().getMap());
+                cf.addMap(dss, rsq, xref.getMap().getMap(), xref.getMap()
+                        .getMappedFromId());
               }
               else
               {
                 // map should be from protein seq to its coding dna
-                cf.addMap(rsq, dss, xref.getMap().getMap().getInverse());
+                cf.addMap(rsq, dss, xref.getMap().getMap().getInverse(),
+                        xref.getMap().getMappedFromId());
               }
             }
           }
@@ -319,153 +333,161 @@ public class CrossRef
        */
       if (!sourceRefs.isEmpty())
       {
-        ASequenceFetcher sftch = SequenceFetcherFactory
-                .getSequenceFetcher();
-        SequenceI[] retrieved = null;
-        try
-        {
-          retrieved = sftch.getSequences(sourceRefs, !fromDna);
-        } catch (Exception e)
-        {
-          System.err
-                  .println("Problem whilst retrieving cross references for Sequence : "
-                          + seq.getName());
-          e.printStackTrace();
-        }
+        retrieveCrossRef(sourceRefs, seq, xrfs, fromDna, cf);
+      }
+    }
+
+    Alignment ral = null;
+    if (rseqs.size() > 0)
+    {
+      ral = new Alignment(rseqs.toArray(new SequenceI[rseqs.size()]));
+      if (!cf.isEmpty())
+      {
+        dataset.addCodonFrame(cf);
+      }
+    }
+    return ral;
+  }
+
+  private void retrieveCrossRef(List<DBRefEntry> sourceRefs, SequenceI seq,
+          DBRefEntry[] xrfs, boolean fromDna, AlignedCodonFrame cf)
+  {
+    ASequenceFetcher sftch = SequenceFetcherFactory.getSequenceFetcher();
+    SequenceI[] retrieved = null;
+    SequenceI dss = seq.getDatasetSequence() == null ? seq : seq
+            .getDatasetSequence();
+    try
+    {
+      retrieved = sftch.getSequences(sourceRefs, !fromDna);
+    } catch (Exception e)
+    {
+      System.err
+              .println("Problem whilst retrieving cross references for Sequence : "
+                      + seq.getName());
+      e.printStackTrace();
+    }
 
-        if (retrieved != null)
+    if (retrieved != null)
+    {
+      updateDbrefMappings(seq, xrfs, retrieved, cf, fromDna);
+      for (SequenceI retrievedSequence : retrieved)
+      {
+        // dataset gets contaminated ccwith non-ds sequences. why ??!
+        // try: Ensembl -> Nuc->Ensembl, Nuc->Uniprot-->Protein->EMBL->
+        SequenceI retrievedDss = retrievedSequence.getDatasetSequence() == null ? retrievedSequence
+                : retrievedSequence.getDatasetSequence();
+        DBRefEntry[] dbr = retrievedSequence.getDBRefs();
+        if (dbr != null)
         {
-          updateDbrefMappings(seq, xrfs, retrieved, cf, fromDna);
-          for (SequenceI retrievedSequence : retrieved)
+          for (DBRefEntry dbref : dbr)
           {
-            SequenceI retrievedDss = retrievedSequence.getDatasetSequence() == null ? retrievedSequence
-                    : retrievedSequence.getDatasetSequence();
-            DBRefEntry[] dbr = retrievedSequence.getDBRefs();
-            if (dbr != null)
+            // find any entry where we should put in the sequence being
+            // cross-referenced into the map
+            Mapping map = dbref.getMap();
+            if (map != null)
             {
-              for (DBRefEntry dbref : dbr)
+              if (map.getTo() != null && map.getMap() != null)
               {
-                // find any entry where we should put in the sequence being
-                // cross-referenced into the map
-                Mapping map = dbref.getMap();
-                if (map != null)
+                // TODO findInDataset requires exact sequence match but
+                // 'congruent' test is only for the mapped part
+                // maybe not a problem in practice since only ENA provide a
+                // mapping and it is to the full protein translation of CDS
+                SequenceI matched = findInDataset(dbref);
+                // matcher.findIdMatch(map.getTo());
+                if (matched != null)
                 {
-                  if (map.getTo() != null && map.getMap() != null)
+                  /*
+                   * already got an xref to this sequence; update this
+                   * map to point to the same sequence, and add
+                   * any new dbrefs to it
+                   */
+                  DBRefEntry[] toRefs = map.getTo().getDBRefs();
+                  if (toRefs != null)
                   {
-                    // TODO findInDataset requires exact sequence match but
-                    // 'congruent' test is only for the mapped part
-                    // maybe not a problem in practice since only ENA provide a
-                    // mapping and it is to the full protein translation of CDS
-                    SequenceI matched = findInDataset(dbref);
-                    // matcher.findIdMatch(map.getTo());
-                    if (matched != null)
-                    {
-                      /*
-                       * already got an xref to this sequence; update this
-                       * map to point to the same sequence, and add
-                       * any new dbrefs to it
-                       */
-                      DBRefEntry[] toRefs = map.getTo().getDBRefs();
-                      if (toRefs != null)
-                      {
-                        for (DBRefEntry ref : toRefs)
-                        {
-                          matched.addDBRef(ref); // add or update mapping
-                        }
-                      }
-                      map.setTo(matched);
-                    }
-                    else
+                    for (DBRefEntry ref : toRefs)
                     {
-                      matcher.add(map.getTo());
+                      matched.addDBRef(ref); // add or update mapping
                     }
-                    try
+                  }
+                  map.setTo(matched);
+                }
+                else
+                {
+                  matcher.add(map.getTo());
+                }
+                try
+                {
+                  // compare ms with dss and replace with dss in mapping
+                  // if map is congruent
+                  SequenceI ms = map.getTo();
+                  int sf = map.getMap().getToLowest();
+                  int st = map.getMap().getToHighest();
+                  SequenceI mappedrg = ms.getSubSequence(sf, st);
+                  // SequenceI loc = dss.getSubSequence(sf, st);
+                  if (mappedrg.getLength() > 0
+                          && ms.getSequenceAsString().equals(
+                                  dss.getSequenceAsString()))
+                  // && mappedrg.getSequenceAsString().equals(
+                  // loc.getSequenceAsString()))
+                  {
+                    String msg = "Mapping updated from " + ms.getName()
+                            + " to retrieved crossreference "
+                            + dss.getName();
+                    System.out.println(msg);
+                    map.setTo(dss);
+
+                    /*
+                     * give the reverse reference the inverse mapping 
+                     * (if it doesn't have one already)
+                     */
+                    setReverseMapping(dss, dbref, cf);
+
+                    /*
+                     * copy sequence features as well, avoiding
+                     * duplication (e.g. same variation from two 
+                     * transcripts)
+                     */
+                    SequenceFeature[] sfs = ms.getSequenceFeatures();
+                    if (sfs != null)
                     {
-                      // compare ms with dss and replace with dss in mapping
-                      // if map is congruent
-                      SequenceI ms = map.getTo();
-                      int sf = map.getMap().getToLowest();
-                      int st = map.getMap().getToHighest();
-                      SequenceI mappedrg = ms.getSubSequence(sf, st);
-                      // SequenceI loc = dss.getSubSequence(sf, st);
-                      if (mappedrg.getLength() > 0
-                              && ms.getSequenceAsString().equals(
-                                      dss.getSequenceAsString()))
-                      // && mappedrg.getSequenceAsString().equals(
-                      // loc.getSequenceAsString()))
+                      for (SequenceFeature feat : sfs)
                       {
-                        String msg = "Mapping updated from " + ms.getName()
-                                + " to retrieved crossreference "
-                                + dss.getName();
-                        System.out.println(msg);
-                        map.setTo(dss);
-
                         /*
-                         * give the reverse reference the inverse mapping 
-                         * (if it doesn't have one already)
+                         * make a flyweight feature object which ignores Parent
+                         * attribute in equality test; this avoids creating many
+                         * otherwise duplicate exon features on genomic sequence
                          */
-                        setReverseMapping(dss, dbref, cf);
-
-                        /*
-                         * copy sequence features as well, avoiding
-                         * duplication (e.g. same variation from two 
-                         * transcripts)
-                         */
-                        SequenceFeature[] sfs = ms.getSequenceFeatures();
-                        if (sfs != null)
+                        SequenceFeature newFeature = new SequenceFeature(
+                                feat)
                         {
-                          for (SequenceFeature feat : sfs)
+                          @Override
+                          public boolean equals(Object o)
                           {
-                            /*
-                             * make a flyweight feature object which ignores Parent
-                             * attribute in equality test; this avoids creating many
-                             * otherwise duplicate exon features on genomic sequence
-                             */
-                            SequenceFeature newFeature = new SequenceFeature(
-                                    feat)
-                            {
-                              @Override
-                              public boolean equals(Object o)
-                              {
-                                return super.equals(o, true);
-                              }
-                            };
-                            dss.addSequenceFeature(newFeature);
+                            return super.equals(o, true);
                           }
-                        }
+                        };
+                        dss.addSequenceFeature(newFeature);
                       }
-                      cf.addMap(retrievedDss, map.getTo(), map.getMap());
-                    } catch (Exception e)
-                    {
-                      System.err
-                              .println("Exception when consolidating Mapped sequence set...");
-                      e.printStackTrace(System.err);
                     }
                   }
+                  cf.addMap(retrievedDss, map.getTo(), map.getMap());
+                } catch (Exception e)
+                {
+                  System.err
+                          .println("Exception when consolidating Mapped sequence set...");
+                  e.printStackTrace(System.err);
                 }
               }
             }
-            retrievedSequence.updatePDBIds();
-            rseqs.add(retrievedDss);
-            dataset.addSequence(retrievedDss);
-            matcher.add(retrievedDss);
           }
         }
+        retrievedSequence.updatePDBIds();
+        rseqs.add(retrievedDss);
+        dataset.addSequence(retrievedDss);
+        matcher.add(retrievedDss);
       }
     }
-
-    Alignment ral = null;
-    if (rseqs.size() > 0)
-    {
-      ral = new Alignment(rseqs.toArray(new SequenceI[rseqs.size()]));
-      if (!cf.isEmpty())
-      {
-        dataset.addCodonFrame(cf);
-      }
-    }
-    return ral;
   }
-
   /**
    * Sets the inverse sequence mapping in the corresponding dbref of the mapped
    * to sequence (if any). This is used after fetching a cross-referenced
@@ -596,14 +618,14 @@ public class CrossRef
   void updateDbrefMappings(SequenceI mapFrom, DBRefEntry[] xrefs,
           SequenceI[] retrieved, AlignedCodonFrame acf, boolean fromDna)
   {
-    SequenceIdMatcher matcher = new SequenceIdMatcher(retrieved);
+    SequenceIdMatcher idMatcher = new SequenceIdMatcher(retrieved);
     for (DBRefEntry xref : xrefs)
     {
       if (!xref.hasMap())
       {
         String targetSeqName = xref.getSource() + "|"
                 + xref.getAccessionId();
-        SequenceI[] matches = matcher.findAllIdMatches(targetSeqName);
+        SequenceI[] matches = idMatcher.findAllIdMatches(targetSeqName);
         if (matches == null)
         {
           return;
@@ -680,6 +702,20 @@ public class CrossRef
       return false;
     }
     xref.setMap(new Mapping(mapTo, mapping));
+
+    /*
+     * and add a reverse DbRef with the inverse mapping
+     */
+    if (mapFrom.getDatasetSequence() != null
+            && mapFrom.getDatasetSequence().getSourceDBRef() != null)
+    {
+      DBRefEntry dbref = new DBRefEntry(mapFrom.getDatasetSequence()
+              .getSourceDBRef());
+      dbref.setMap(new Mapping(mapFrom.getDatasetSequence(), mapping
+              .getInverse()));
+      mapTo.addDBRef(dbref);
+    }
+
     if (fromDna)
     {
       AlignmentUtils.computeProteinFeatures(mapFrom, mapTo, mapping);
@@ -703,11 +739,11 @@ public class CrossRef
    *          context was searching from Protein sequences
    * @param sequenceI
    * @param lrfs
-   * @param rseqs
+   * @param foundSeqs
    * @return true if matches were found.
    */
   private boolean searchDatasetXrefs(boolean fromDna, SequenceI sequenceI,
-          DBRefEntry[] lrfs, List<SequenceI> rseqs, AlignedCodonFrame cf)
+          DBRefEntry[] lrfs, List<SequenceI> foundSeqs, AlignedCodonFrame cf)
   {
     boolean found = false;
     if (lrfs == null)
@@ -720,7 +756,7 @@ public class CrossRef
       // add in wildcards
       xref.setVersion(null);
       xref.setMap(null);
-      found |= searchDataset(fromDna, sequenceI, xref, rseqs, cf, false);
+      found |= searchDataset(fromDna, sequenceI, xref, foundSeqs, cf, false);
     }
     return found;
   }
@@ -732,20 +768,29 @@ public class CrossRef
    * @param fromDna
    *          true if context was searching for refs *from* dna sequence, false
    *          if context was searching for refs *from* protein sequence
-   * @param sequenceI
+   * @param fromSeq
    *          a sequence to ignore (start point of search)
    * @param xrf
    *          a cross-reference to try to match
-   * @param rseqs
+   * @param foundSeqs
    *          result list to add to
-   * @param cf
+   * @param mappings
    *          a set of sequence mappings to add to
    * @param direct
-   *          - search all references or only subset
+   *          - indicates the type of relationship between returned sequences,
+   *          xrf, and sequenceI that is required.
+   *          <ul>
+   *          <li>direct implies xrf is a primary reference for sequenceI AND
+   *          the sequences to be located (eg a uniprot ID for a protein
+   *          sequence, and a uniprot ref on a transcript sequence).</li>
+   *          <li>indirect means xrf is a cross reference with respect to
+   *          sequenceI or all the returned sequences (eg a genomic reference
+   *          associated with a locus and one or more transcripts)</li>
+   *          </ul>
    * @return true if relationship found and sequence added.
    */
-  boolean searchDataset(boolean fromDna, SequenceI sequenceI,
-          DBRefEntry xrf, List<SequenceI> rseqs, AlignedCodonFrame cf,
+  boolean searchDataset(boolean fromDna, SequenceI fromSeq,
+          DBRefEntry xrf, List<SequenceI> foundSeqs, AlignedCodonFrame mappings,
           boolean direct)
   {
     boolean found = false;
@@ -768,9 +813,13 @@ public class CrossRef
           if (nxt.getDatasetSequence() != null)
           {
             System.err
-                    .println("Implementation warning: getProducts passed a dataset alignment without dataset sequences in it!");
+                    .println("Implementation warning: CrossRef initialised with a dataset alignment with non-dataset sequences in it! ("
+                            + nxt.getDisplayId(true)
+                            + " has ds reference "
+                            + nxt.getDatasetSequence().getDisplayId(true)
+                            + ")");
           }
-          if (nxt == sequenceI || nxt == sequenceI.getDatasetSequence())
+          if (nxt == fromSeq || nxt == fromSeq.getDatasetSequence())
           {
             continue;
           }
@@ -779,8 +828,7 @@ public class CrossRef
            * complementary type if !direct
            */
           {
-            boolean isDna = Comparison
-                    .isNucleotide(new SequenceI[] { nxt });
+            boolean isDna = !nxt.isProtein();
             if (direct ? (isDna != fromDna) : (isDna == fromDna))
             {
               // skip this sequence because it is wrong molecule type
@@ -791,30 +839,29 @@ public class CrossRef
           // look for direct or indirect references in common
           DBRefEntry[] poss = nxt.getDBRefs();
           List<DBRefEntry> cands = null;
-          /*
-           * TODO does this make any sense?
-           * if 'direct', search the dbrefs for xrf
-           * else, filter the dbrefs by type and then search for xrf
-           * - the result is the same isn't it?
-           */
-          if (direct)
-          {
-            cands = DBRefUtils.searchRefs(poss, xrf);
-          }
-          else
-          {
-            poss = DBRefUtils.selectDbRefs(!fromDna, poss);
-            cands = DBRefUtils.searchRefs(poss, xrf);
-          }
+
+          // todo: indirect specifies we select either direct references to nxt
+          // that match xrf which is indirect to sequenceI, or indirect
+          // references to nxt that match xrf which is direct to sequenceI
+          cands = DBRefUtils.searchRefs(poss, xrf);
+          // else
+          // {
+          // poss = DBRefUtils.selectDbRefs(nxt.isProtein()!fromDna, poss);
+          // cands = DBRefUtils.searchRefs(poss, xrf);
+          // }
           if (!cands.isEmpty())
           {
-            if (!rseqs.contains(nxt))
+            if (!foundSeqs.contains(nxt))
             {
               found = true;
-              rseqs.add(nxt);
-              if (cf != null)
+              foundSeqs.add(nxt);
+              if (mappings != null && !direct)
               {
-                // don't search if we aren't given a codon map object
+                /*
+                 * if the matched sequence has mapped dbrefs to
+                 * protein product / cdna, add equivalent mappings to
+                 * our source sequence
+                 */
                 for (DBRefEntry candidate : cands)
                 {
                   Mapping mapping = candidate.getMap();
@@ -824,23 +871,21 @@ public class CrossRef
                     if (mapping.getTo() != null
                             && map.getFromRatio() != map.getToRatio())
                     {
-                      // get sense of map correct for adding to product
-                      // alignment.
-                      if (fromDna)
+                      /*
+                       * add a mapping, as from dna to peptide sequence
+                       */
+                      if (map.getFromRatio() == 3)
                       {
-                        // map is from dna seq to a protein product
-                        cf.addMap(sequenceI, nxt, map);
+                        mappings.addMap(nxt, fromSeq, map);
                       }
                       else
                       {
-                        // map should be from protein seq to its coding dna
-                        cf.addMap(nxt, sequenceI, map.getInverse());
+                        mappings.addMap(nxt, fromSeq, map.getInverse());
                       }
                     }
                   }
                 }
               }
-              // TODO: add mapping between sequences if necessary
             }
           }
         }