JAL-2781 mark aligned positions via StructureMapping’s Mapping object, and transfer...
authorJim Procter <jprocter@issues.jalview.org>
Sun, 11 Feb 2018 10:54:00 +0000 (10:54 +0000)
committerJim Procter <jprocter@issues.jalview.org>
Sun, 11 Feb 2018 10:59:14 +0000 (10:59 +0000)
src/MCview/PDBChain.java
src/jalview/structure/StructureMapping.java
src/jalview/ws/sifts/SiftsClient.java

index f4bd31c..29b994e 100755 (executable)
@@ -162,6 +162,50 @@ public class PDBChain
   }
 
   /**
+   * Annotate the residues with their corresponding positions in s1 using the
+   * alignment in as NOTE: This clears all atom.alignmentMapping values on the
+   * structure.
+   * 
+   * @param as
+   * @param s1
+   */
+  public void makeExactMapping(StructureMapping mapping, SequenceI s1)
+  {
+    // first clear out any old alignmentMapping values:
+    for (Atom atom : atoms)
+    {
+      atom.alignmentMapping = -1;
+    }
+    SequenceI ds = s1;
+    while (ds.getDatasetSequence() != null)
+    {
+      ds = ds.getDatasetSequence();
+    }
+    int pdboffset = 0;
+    for (Residue res : residues)
+    {
+      // res.number isn't set correctly for discontinuous/mismapped residues
+      int seqpos = mapping.getSeqPos(res.atoms.get(0).resNumber);
+      char strchar = sequence.getCharAt(pdboffset++);
+      if (seqpos == StructureMapping.UNASSIGNED_VALUE)
+      {
+        continue;
+      }
+      char seqchar = ds.getCharAt(seqpos - ds.getStart());
+
+      boolean sameResidue = Comparison.isSameResidue(
+              seqchar, strchar, false);
+      if (sameResidue)
+      {
+        for (Atom atom : res.atoms)
+        {
+          atom.alignmentMapping = seqpos - 1;
+        }
+      }
+    }
+  }
+
+  /**
    * Copies over the RESNUM seqfeatures from the internal chain sequence to the
    * mapped sequence
    * 
@@ -550,6 +594,12 @@ public class PDBChain
   {
     SequenceI sq = mapping.getSequence();
     SequenceI dsq = sq;
+    if (sqmpping == null)
+    {
+      // SIFTS mappings are recorded in the StructureMapping object...
+
+      sqmpping = mapping.getSeqToPdbMapping();
+    }
     if (sq != null)
     {
       while (dsq.getDatasetSequence() != null)
index c103521..fcf322d 100644 (file)
@@ -21,6 +21,7 @@
 package jalview.structure;
 
 import jalview.datamodel.AlignmentAnnotation;
+import jalview.datamodel.Mapping;
 import jalview.datamodel.SequenceI;
 
 import java.util.ArrayList;
@@ -49,6 +50,7 @@ public class StructureMapping
   // and atomNo
   HashMap<Integer, int[]> mapping;
 
+  jalview.datamodel.Mapping seqToPdbMapping = null;
   /**
    * Constructor
    * 
@@ -73,6 +75,14 @@ public class StructureMapping
     this.mappingDetails = mappingDetails;
   }
 
+  public StructureMapping(SequenceI seq, String pdbFile2, String pdbId2,
+          String chain, HashMap<Integer, int[]> mapping2,
+          String mappingOutput, Mapping seqToPdbMapping)
+  {
+    this(seq, pdbFile2, pdbId2, chain, mapping2, mappingOutput);
+    this.seqToPdbMapping = seqToPdbMapping;
+  }
+
   public SequenceI getSequence()
   {
     return sequence;
@@ -248,4 +258,9 @@ public class StructureMapping
   {
     return mapping;
   }
+
+  public Mapping getSeqToPdbMapping()
+  {
+    return seqToPdbMapping;
+  }
 }
index 289030e..e1b539e 100644 (file)
@@ -448,8 +448,6 @@ public class SiftsClient implements SiftsClientI
     HashMap<Integer, int[]> mapping = getGreedyMapping(chain, seq, ps);
 
     String mappingOutput = mappingDetails.toString();
-    StructureMapping siftsMapping = new StructureMapping(seq, pdbFile,
-            pdbId, chain, mapping, mappingOutput);
     StructureMapping siftsMapping = new StructureMapping(aseq, pdbFile,
             pdbId, chain, mapping, mappingOutput, seqFromPdbMapping);
 
@@ -514,15 +512,61 @@ public class SiftsClient implements SiftsClientI
     {
       throw new SiftsException("SIFTS mapping failed");
     }
+    // also construct a mapping object between the seq-coord sys and the PDB seq's coord sys
 
     Integer[] keys = mapping.keySet().toArray(new Integer[0]);
     Arrays.sort(keys);
     seqStart = keys[0];
     seqEnd = keys[keys.length - 1];
-
+    List<int[]> from=new ArrayList<>(),to=new ArrayList<>();
+    int[]_cfrom=null,_cto=null;
     String matchedSeq = originalSeq;
-    if (seqStart != UNASSIGNED)
+    if (seqStart != UNASSIGNED) // fixme! seqStart can map to -1 for a pdb sequence that starts <-1
     {
+      for (int seqps:keys)
+      {
+        int pdbpos = mapping.get(seqps)[PDBE_POS];
+        if (pdbpos == UNASSIGNED)
+        {
+          // not correct - pdbpos might be -1, but leave it for now
+          continue;
+        }
+        if (_cfrom==null || seqps!=_cfrom[1]+1)
+        {
+          _cfrom = new int[] { seqps,seqps};
+          from.add(_cfrom);
+          _cto = null; // discontinuity
+        } else {
+          _cfrom[1]= seqps;
+        }
+        if (_cto==null || pdbpos!=1+_cto[1])
+        {
+          _cto = new int[] { pdbpos,pdbpos};
+          to.add(_cto);
+        } else {
+          _cto[1] = pdbpos;
+        }
+      }
+      _cfrom = new int[from.size() * 2];
+      _cto = new int[to.size() * 2];
+      int p = 0;
+      for (int[] range : from)
+      {
+        _cfrom[p++] = range[0];
+        _cfrom[p++] = range[1];
+      }
+      ;
+      p = 0;
+      for (int[] range : to)
+      {
+        _cto[p++] = range[0];
+        _cto[p++] = range[1];
+      }
+      ;
+
+      seqFromPdbMapping = new jalview.datamodel.Mapping(null, _cto, _cfrom,
+              1,
+              1);
       pdbStart = mapping.get(seqStart)[PDB_RES_POS];
       pdbEnd = mapping.get(seqEnd)[PDB_RES_POS];
       int orignalSeqStart = seq.getStart();
@@ -575,6 +619,8 @@ public class SiftsClient implements SiftsClientI
     TreeMap<Integer, String> resNumMap = shp.getResNumMap();
     List<Integer> omitNonObserved = shp.getOmitNonObserved();
     int nonObservedShiftIndex = shp.getNonObservedShiftIndex();
+    int pdbeNonObservedCount = shp.getPdbeNonObserved();
+    int firstPDBResNum = UNASSIGNED;
     for (Segment segment : segments)
     {
       // System.out.println("Mapping segments : " + segment.getSegId() + "\\"s