src/jalview/ws/datamodel/alphafold/MappableContactMatrix.java

   1 package jalview.ws.datamodel.alphafold;
   2
   3 import java.util.ArrayList;
   4
   5 import jalview.datamodel.ContactListI;
   6 import jalview.datamodel.ContactListImpl;
   7 import jalview.datamodel.ContactListProviderI;
   8 import jalview.datamodel.Mapping;
   9 import jalview.datamodel.SequenceI;
  10 import jalview.util.MapList;
  11 import jalview.ws.datamodel.MappableContactMatrixI;
  12
  13 public abstract class MappableContactMatrix<T extends MappableContactMatrix<T>> implements MappableContactMatrixI
  14 {
  15   SequenceI refSeq = null;
  16   MapList toSeq = null;
  17
  18   /**
  19    * the length that refSeq is expected to be (excluding gaps, of course)
  20    */
  21   int length;
  22
  23
  24   @Override
  25   public boolean hasReferenceSeq()
  26   {
  27     return (refSeq != null);
  28   }
  29
  30   @Override
  31   public SequenceI getReferenceSeq()
  32   {
  33     return refSeq;
  34   }
  35   @Override
  36   public MapList getMapFor(SequenceI mapSeq)
  37   {
  38     if (refSeq!=null)
  39     {
  40       while (mapSeq!=refSeq && mapSeq.getDatasetSequence()!=null)
  41       {
  42         mapSeq = mapSeq.getDatasetSequence();
  43       }
  44       if (mapSeq!=refSeq)
  45       {
  46         return null;
  47       }
  48     } else {
  49       if (mapSeq!=null) {
  50         // our MapList does not concern this seq
  51         return null;
  52       }
  53     }
  54
  55     return toSeq;
  56   }
  57
  58   protected void setRefSeq(SequenceI _refSeq)
  59   {
  60     refSeq = _refSeq;
  61     while (refSeq.getDatasetSequence() != null)
  62     {
  63       refSeq = refSeq.getDatasetSequence();
  64     }
  65     length = _refSeq.getEnd() - _refSeq.getStart() + 1;
  66 //    if (length!=refSeq.getLength() || _refSeq.getStart()!=1)
  67     {
  68       toSeq = new MapList(new int[] { _refSeq.getStart(), _refSeq.getEnd()}, new int[] { 0,length-1}, 1,1);
  69     }
  70   }
  71
  72   public T liftOver(SequenceI newRefSeq, Mapping sp2sq)
  73   {
  74     if (sp2sq.getMappedWidth() != sp2sq.getWidth())
  75     {
  76       // TODO: employ getWord/MappedWord to transfer annotation between cDNA and
  77       // Protein reference frames
  78       throw new Error(
  79               "liftOver currently not implemented for transfer of annotation between different types of seqeunce");
  80     }
  81     boolean mapIsTo = (sp2sq != null) ? (sp2sq.getTo() == refSeq) : false;
  82
  83     /**
  84      * map from matrix to toSeq's coordinate frame
  85      */
  86     int[] refMap = toSeq.locateInFrom(0, length - 1);
  87     ArrayList<Integer> newFromMap = new ArrayList<Integer>();
  88     int last = -1;
  89     for (int i = 0; i < refMap.length; i += 2)
  90     {
  91       /*
  92        * for each contiguous range in toSeq, locate corresponding range in sequence mapped to toSeq by sp2sq
  93        */
  94       int[] sp2map = mapIsTo
  95               ? sp2sq.getMap().locateInFrom(refMap[i], refMap[i + 1])
  96               : sp2sq.getMap().locateInTo(refMap[i], refMap[i + 1]);
  97       if (sp2map == null)
  98       {
  99         continue;
 100       }
 101
 102       for (int spm = 0; spm < sp2map.length; spm += 2)
 103       {
 104
 105         if (last > -1)
 106         {
 107           if (sp2map[spm] != last + 1)
 108           {
 109             newFromMap.add(sp2map[spm]);
 110           }
 111           else
 112           {
 113             newFromMap.remove(newFromMap.size() - 1);
 114           }
 115         }
 116         else
 117         {
 118           newFromMap.add(sp2map[spm]);
 119         }
 120         last = sp2map[spm + 1];
 121         newFromMap.add(last);
 122       }
 123     }
 124     if ((newFromMap.size() % 2) != 0)
 125     {
 126       // should have had an even number of int ranges!
 127       throw new Error("PAEMatrix liftover failed.");
 128     }
 129     int fromIntMap[] = new int[newFromMap.size()];
 130     int ipos = 0;
 131     for (Integer i : newFromMap)
 132     {
 133       fromIntMap[ipos++] = i;
 134     }
 135     MapList newFromMapList = new MapList(fromIntMap,
 136             new int[]
 137             { 0, length - 1 }, 1, 1);
 138
 139     T newCM = newMappableContactMatrix(newRefSeq, newFromMapList);
 140     return newCM;
 141   }
 142
 143   protected abstract T  newMappableContactMatrix(SequenceI newRefSeq,
 144           MapList newFromMapList);
 145   @Override
 146   public ContactListI getMappableContactList(final SequenceI localFrame,
 147           final int column)
 148   {
 149     final int _column;
 150     final int _lcolumn;
 151     if (localFrame==null)
 152     {
 153       throw new Error("Unimplemented when no local sequence given.");
 154     }
 155     // return a ContactListI for column
 156     // column is index into localFrame
 157     // 1. map column to corresponding column in matrix
 158
 159     _lcolumn=localFrame.findPosition(column);
 160
 161     if (toSeq != null)
 162     {
 163       SequenceI lf = localFrame, uf = refSeq;
 164
 165       // just look for dataset sequences and check they are the same.
 166       // in future we could use DBRefMappings/whatever.
 167       while (lf.getDatasetSequence() != null
 168               || uf.getDatasetSequence() != null)
 169       {
 170         if (lf.getDatasetSequence() != null)
 171         {
 172           lf = lf.getDatasetSequence();
 173         }
 174         if (uf.getDatasetSequence() != null)
 175         {
 176           uf = uf.getDatasetSequence();
 177         }
 178       }
 179       if (lf != uf)
 180       {
 181         // could try harder to find a mapping
 182         throw new Error("This Matrix associated with '" + refSeq.getName()
 183                 + "' is not mappable for the given localFrame sequence. ("
 184                 + localFrame.getName() + ")");
 185       }
 186       // check the mapping to see if localFrame _lcolumn exists
 187       int[] word = toSeq.locateInTo(_lcolumn, _lcolumn);
 188       if (word == null)
 189       {
 190         return null;
 191       }
 192       _column = word[0];
 193     }
 194     else
 195     {
 196       // no mapping
 197       _column = _lcolumn;
 198     }
 199
 200     // TODO - remove ? this may be a redundant check
 201     if (_column < 0 || ((toSeq != null && _column > toSeq.getToHighest())
 202             || (toSeq == null && getHeight() <= _column)))
 203     {
 204       return null;
 205     }
 206
 207     // 2. resolve ranges in matrix corresponding to range in localFrame
 208     final int[] matrixRange = toSeq == null
 209             ? new int[]
 210             { localFrame.getStart(), localFrame.getEnd() }
 211             : toSeq.locateInTo(localFrame.getStart(), localFrame.getEnd());
 212
 213     int h = 0;
 214     for (int p = 0; p < matrixRange.length; p += 2)
 215     {
 216       h += 1+Math.abs(matrixRange[p + 1] - matrixRange[p]);
 217     }
 218     final int rangeHeight = h;
 219     // 3. Construct ContactListImpl instance for just those segments.
 220
 221     return new ContactListImpl(new ContactListProviderI()
 222     {
 223
 224       public int getColumn()
 225       {
 226         return column;
 227       }
 228       @Override
 229       public int getPosition()
 230       {
 231         return _column;
 232       }
 233
 234       @Override
 235       public int getContactHeight()
 236       {
 237         return rangeHeight;
 238       }
 239
 240       @Override
 241       public double getContactAt(int mcolumn)
 242       {
 243         if (mcolumn<0 || mcolumn>=rangeHeight)
 244         {
 245           return -1;
 246         }
 247         return getElementAt(_column, locateInRange(mcolumn));
 248
 249         // this code maps from mcolumn to localFrame - but that isn't what's needed
 250 //        int loccolumn = localFrame.findPosition(mcolumn);
 251 //        int[] lcolumn=(toSeq==null) ? new int[] {mcolumn} : toSeq.locateInTo(loccolumn,loccolumn);
 252 //        if (lcolumn==null || lcolumn[0] < 0 || lcolumn[0] >= rangeHeight)
 253 //        {
 254 //          return -1;
 255 //        }
 256 //        return getElementAt(_column,lcolumn[0]);
 257       }
 258       /**
 259        * @return the mcolumn'th position in the matrixRange window on the matrix
 260        */
 261       private int locateInRange(int mcolumn)
 262       {
 263
 264         int h=0,p=0;
 265         while (h < mcolumn && p+2 < matrixRange.length)
 266         {
 267           h += 1+Math.abs(matrixRange[p + 1] - matrixRange[p]);
 268           p+=2;
 269         }
 270         return matrixRange[p]+mcolumn-h;
 271       }
 272
 273       @Override
 274       public int[] getMappedPositionsFor(int cStart, int cEnd)
 275       {
 276         if (!hasReferenceSeq())
 277         {
 278           return ContactListProviderI.super.getMappedPositionsFor(cStart, cEnd);
 279         }
 280         // map into segment of matrix being shown
 281         int realCstart = locateInRange(cStart);
 282         int realCend = locateInRange(cEnd);
 283
 284         // TODO account for discontinuities in the mapping
 285
 286         int[] mappedPositions = toSeq.locateInFrom(realCstart,realCend);
 287         if (mappedPositions!=null) {
 288           int s=-1,e=-1;
 289           for (int p=0;p<mappedPositions.length;p++)
 290           {
 291             if (s==-1 && mappedPositions[p]>=localFrame.getStart())
 292             {
 293               s=p; // remember first position within local frame
 294             }
 295             if (e==-1 || mappedPositions[p]<=localFrame.getEnd())
 296             {
 297               // update end pointer
 298               e=p;
 299               // compute local map
 300               mappedPositions[p] = localFrame.findIndex(mappedPositions[p]);
 301             }
 302           }
 303         }
 304         return mappedPositions;
 305       }
 306     });
 307   }
 308
 309   /**
 310    * get a specific element of the contact matrix in its data-local coordinates
 311    * rather than the mapped frame. Implementations are allowed to throw RunTimeExceptions if _column/i are out of bounds
 312    *
 313    * @param _column
 314    * @param i
 315    * @return
 316    */
 317   protected abstract double getElementAt(int _column, int i);
 318 }