src/jalview/ws/sifts/SiftsClient.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.ws.sifts;
  22
  23 import java.util.Locale;
  24
  25 import java.io.File;
  26 import java.io.FileInputStream;
  27 import java.io.FileOutputStream;
  28 import java.io.IOException;
  29 import java.io.InputStream;
  30 import java.io.PrintStream;
  31 import java.net.URL;
  32 import java.net.URLConnection;
  33 import java.nio.file.Files;
  34 import java.nio.file.Path;
  35 import java.nio.file.attribute.BasicFileAttributes;
  36 import java.util.ArrayList;
  37 import java.util.Arrays;
  38 import java.util.Collection;
  39 import java.util.Collections;
  40 import java.util.Date;
  41 import java.util.HashMap;
  42 import java.util.HashSet;
  43 import java.util.List;
  44 import java.util.Map;
  45 import java.util.Set;
  46 import java.util.TreeMap;
  47 import java.util.zip.GZIPInputStream;
  48
  49 import javax.xml.bind.JAXBContext;
  50 import javax.xml.bind.JAXBElement;
  51 import javax.xml.bind.Unmarshaller;
  52 import javax.xml.stream.XMLInputFactory;
  53 import javax.xml.stream.XMLStreamReader;
  54
  55 import jalview.analysis.AlignSeq;
  56 import jalview.analysis.scoremodels.ScoreMatrix;
  57 import jalview.analysis.scoremodels.ScoreModels;
  58 import jalview.api.DBRefEntryI;
  59 import jalview.api.SiftsClientI;
  60 import jalview.bin.Console;
  61 import jalview.datamodel.DBRefEntry;
  62 import jalview.datamodel.DBRefSource;
  63 import jalview.datamodel.SequenceI;
  64 import jalview.io.BackupFiles;
  65 import jalview.io.StructureFile;
  66 import jalview.schemes.ResidueProperties;
  67 import jalview.structure.StructureMapping;
  68 import jalview.util.Comparison;
  69 import jalview.util.DBRefUtils;
  70 import jalview.util.Format;
  71 import jalview.util.Platform;
  72 import jalview.xml.binding.sifts.Entry;
  73 import jalview.xml.binding.sifts.Entry.Entity;
  74 import jalview.xml.binding.sifts.Entry.Entity.Segment;
  75 import jalview.xml.binding.sifts.Entry.Entity.Segment.ListMapRegion.MapRegion;
  76 import jalview.xml.binding.sifts.Entry.Entity.Segment.ListResidue.Residue;
  77 import jalview.xml.binding.sifts.Entry.Entity.Segment.ListResidue.Residue.CrossRefDb;
  78 import jalview.xml.binding.sifts.Entry.Entity.Segment.ListResidue.Residue.ResidueDetail;
  79 import mc_view.Atom;
  80 import mc_view.PDBChain;
  81
  82 public class SiftsClient implements SiftsClientI
  83 {
  84   /*
  85    * for use in mocking out file fetch for tests only
  86    * - reset to null after testing!
  87    */
  88   private static File mockSiftsFile;
  89
  90   private Entry siftsEntry;
  91
  92   private StructureFile pdb;
  93
  94   private String pdbId;
  95
  96   private String structId;
  97
  98   private CoordinateSys seqCoordSys = CoordinateSys.UNIPROT;
  99
 100   /**
 101    * PDB sequence position to sequence coordinate mapping as derived from SIFTS
 102    * record for the identified SeqCoordSys Used for lift-over from sequence
 103    * derived from PDB (with first extracted PDBRESNUM as 'start' to the sequence
 104    * being annotated with PDB data
 105    */
 106   private jalview.datamodel.Mapping seqFromPdbMapping;
 107
 108   private static final int BUFFER_SIZE = 4096;
 109
 110   public static final int UNASSIGNED = Integer.MIN_VALUE;
 111
 112   private static final int PDB_RES_POS = 0;
 113
 114   private static final int PDB_ATOM_POS = 1;
 115
 116   private static final int PDBE_POS = 2;
 117
 118   private static final String NOT_OBSERVED = "Not_Observed";
 119
 120   private static final String SIFTS_SPLIT_FTP_BASE_URL = "https://ftp.ebi.ac.uk/pub/databases/msd/sifts/split_xml/";
 121
 122   private final static String NEWLINE = System.lineSeparator();
 123
 124   private String curSourceDBRef;
 125
 126   private HashSet<String> curDBRefAccessionIdsString;
 127
 128   private enum CoordinateSys
 129   {
 130     UNIPROT("UniProt"), PDB("PDBresnum"), PDBe("PDBe");
 131
 132     private String name;
 133
 134     private CoordinateSys(String name)
 135     {
 136       this.name = name;
 137     }
 138
 139     public String getName()
 140     {
 141       return name;
 142     }
 143   };
 144
 145   private enum ResidueDetailType
 146   {
 147     NAME_SEC_STRUCTURE("nameSecondaryStructure"),
 148     CODE_SEC_STRUCTURE("codeSecondaryStructure"), ANNOTATION("Annotation");
 149
 150     private String code;
 151
 152     private ResidueDetailType(String code)
 153     {
 154       this.code = code;
 155     }
 156
 157     public String getCode()
 158     {
 159       return code;
 160     }
 161   };
 162
 163   /**
 164    * Fetch SIFTs file for the given PDBfile and construct an instance of
 165    * SiftsClient
 166    *
 167    * @param pdbId
 168    * @throws SiftsException
 169    */
 170   public SiftsClient(StructureFile pdb) throws SiftsException
 171   {
 172     this.pdb = pdb;
 173     this.pdbId = pdb.getId();
 174     File siftsFile = getSiftsFile(pdbId);
 175     siftsEntry = parseSIFTs(siftsFile);
 176   }
 177
 178   /**
 179    * Parse the given SIFTs File and return a JAXB POJO of parsed data
 180    *
 181    * @param siftFile
 182    *          - the GZipped SIFTs XML file to parse
 183    * @return
 184    * @throws Exception
 185    *           if a problem occurs while parsing the SIFTs XML
 186    */
 187   private Entry parseSIFTs(File siftFile) throws SiftsException
 188   {
 189     try (InputStream in = new FileInputStream(siftFile);
 190             GZIPInputStream gzis = new GZIPInputStream(in);)
 191     {
 192       // System.out.println("File : " + siftFile.getAbsolutePath());
 193       JAXBContext jc = JAXBContext.newInstance("jalview.xml.binding.sifts");
 194       XMLStreamReader streamReader = XMLInputFactory.newInstance()
 195               .createXMLStreamReader(gzis);
 196       Unmarshaller um = jc.createUnmarshaller();
 197       JAXBElement<Entry> jbe = um.unmarshal(streamReader, Entry.class);
 198       return jbe.getValue();
 199     } catch (Exception e)
 200     {
 201       e.printStackTrace();
 202       throw new SiftsException(e.getMessage());
 203     }
 204   }
 205
 206   /**
 207    * Get a SIFTs XML file for a given PDB Id from Cache or download from FTP
 208    * repository if not found in cache
 209    *
 210    * @param pdbId
 211    * @return SIFTs XML file
 212    * @throws SiftsException
 213    */
 214   public static File getSiftsFile(String pdbId) throws SiftsException
 215   {
 216     /*
 217      * return mocked file if it has been set
 218      */
 219     if (mockSiftsFile != null)
 220     {
 221       return mockSiftsFile;
 222     }
 223
 224     String siftsFileName = SiftsSettings.getSiftDownloadDirectory()
 225             + pdbId.toLowerCase(Locale.ROOT) + ".xml.gz";
 226     File siftsFile = new File(siftsFileName);
 227     if (siftsFile.exists())
 228     {
 229       // The line below is required for unit testing... don't comment it out!!!
 230       System.out.println(">>> SIFTS File already downloaded for " + pdbId);
 231
 232       if (isFileOlderThanThreshold(siftsFile,
 233               SiftsSettings.getCacheThresholdInDays()))
 234       {
 235         File oldSiftsFile = new File(siftsFileName + "_old");
 236         BackupFiles.moveFileToFile(siftsFile, oldSiftsFile);
 237         try
 238         {
 239           siftsFile = downloadSiftsFile(pdbId.toLowerCase(Locale.ROOT));
 240           oldSiftsFile.delete();
 241           return siftsFile;
 242         } catch (IOException e)
 243         {
 244           e.printStackTrace();
 245           BackupFiles.moveFileToFile(oldSiftsFile, siftsFile);
 246           return new File(siftsFileName);
 247         }
 248       }
 249       else
 250       {
 251         return siftsFile;
 252       }
 253     }
 254     try
 255     {
 256       siftsFile = downloadSiftsFile(pdbId.toLowerCase(Locale.ROOT));
 257     } catch (IOException e)
 258     {
 259       throw new SiftsException(e.getMessage());
 260     }
 261     return siftsFile;
 262   }
 263
 264   /**
 265    * This method enables checking if a cached file has exceeded a certain
 266    * threshold(in days)
 267    *
 268    * @param file
 269    *          the cached file
 270    * @param noOfDays
 271    *          the threshold in days
 272    * @return
 273    */
 274   public static boolean isFileOlderThanThreshold(File file, int noOfDays)
 275   {
 276     Path filePath = file.toPath();
 277     BasicFileAttributes attr;
 278     int diffInDays = 0;
 279     try
 280     {
 281       attr = Files.readAttributes(filePath, BasicFileAttributes.class);
 282       diffInDays = (int) ((new Date().getTime()
 283               - attr.lastModifiedTime().toMillis())
 284               / (1000 * 60 * 60 * 24));
 285       // System.out.println("Diff in days : " + diffInDays);
 286     } catch (IOException e)
 287     {
 288       e.printStackTrace();
 289     }
 290     return noOfDays <= diffInDays;
 291   }
 292
 293   /**
 294    * Download a SIFTs XML file for a given PDB Id from an FTP repository
 295    *
 296    * @param pdbId
 297    * @return downloaded SIFTs XML file
 298    * @throws SiftsException
 299    * @throws IOException
 300    */
 301   public static File downloadSiftsFile(String pdbId)
 302           throws SiftsException, IOException
 303   {
 304     if (pdbId.contains(".cif"))
 305     {
 306       pdbId = pdbId.replace(".cif", "");
 307     }
 308     String siftFile = pdbId + ".xml.gz";
 309     String siftsFileFTPURL = getDownloadUrlFor(siftFile);
 310
 311     /*
 312      * Download the file from URL to either
 313      * Java: directory of cached downloaded SIFTS files
 314      * Javascript: temporary 'file' (in-memory cache)
 315      */
 316     File downloadTo = null;
 317     if (Platform.isJS())
 318     {
 319       downloadTo = File.createTempFile(siftFile, ".xml.gz");
 320     }
 321     else
 322     {
 323       downloadTo = new File(
 324               SiftsSettings.getSiftDownloadDirectory() + siftFile);
 325       File siftsDownloadDir = new File(
 326               SiftsSettings.getSiftDownloadDirectory());
 327       if (!siftsDownloadDir.exists())
 328       {
 329         siftsDownloadDir.mkdirs();
 330       }
 331     }
 332
 333     // System.out.println(">> Download ftp url : " + siftsFileFTPURL);
 334     // long now = System.currentTimeMillis();
 335     URL url = new URL(siftsFileFTPURL);
 336     URLConnection conn = url.openConnection();
 337     InputStream inputStream = conn.getInputStream();
 338     FileOutputStream outputStream = new FileOutputStream(downloadTo);
 339     byte[] buffer = new byte[BUFFER_SIZE];
 340     int bytesRead = -1;
 341     while ((bytesRead = inputStream.read(buffer)) != -1)
 342     {
 343       outputStream.write(buffer, 0, bytesRead);
 344     }
 345     outputStream.close();
 346     inputStream.close();
 347     // System.out.println(">>> File downloaded : " + downloadedSiftsFile
 348     // + " took " + (System.currentTimeMillis() - now) + "ms");
 349     return downloadTo;
 350   }
 351
 352   public static String getDownloadUrlFor(String siftFile)
 353   {
 354     String durl = SIFTS_SPLIT_FTP_BASE_URL+siftFile.substring(1, 3)+"/"+siftFile;
 355     Console.trace("SIFTS URL for "+siftFile+" is "+durl);
 356     return durl;
 357
 358   }
 359
 360   /**
 361    * Delete the SIFTs file for the given PDB Id in the local SIFTs download
 362    * directory
 363    *
 364    * @param pdbId
 365    * @return true if the file was deleted or doesn't exist
 366    */
 367   public static boolean deleteSiftsFileByPDBId(String pdbId)
 368   {
 369     File siftsFile = new File(SiftsSettings.getSiftDownloadDirectory()
 370             + pdbId.toLowerCase(Locale.ROOT) + ".xml.gz");
 371     if (siftsFile.exists())
 372     {
 373       return siftsFile.delete();
 374     }
 375     return true;
 376   }
 377
 378   /**
 379    * Get a valid SIFTs DBRef for the given sequence current SIFTs entry
 380    *
 381    * @param seq
 382    *          - the target sequence for the operation
 383    * @return a valid DBRefEntry that is SIFTs compatible
 384    * @throws Exception
 385    *           if no valid source DBRefEntry was found for the given sequences
 386    */
 387   public DBRefEntryI getValidSourceDBRef(SequenceI seq)
 388           throws SiftsException
 389   {
 390     List<DBRefEntry> dbRefs = seq.getPrimaryDBRefs();
 391     if (dbRefs == null || dbRefs.size() < 1)
 392     {
 393       throw new SiftsException(
 394               "Source DBRef could not be determined. DBRefs might not have been retrieved.");
 395     }
 396
 397     for (DBRefEntry dbRef : dbRefs)
 398     {
 399       if (dbRef == null || dbRef.getAccessionId() == null
 400               || dbRef.getSource() == null)
 401       {
 402         continue;
 403       }
 404       String canonicalSource = DBRefUtils
 405               .getCanonicalName(dbRef.getSource());
 406       if (isValidDBRefEntry(dbRef)
 407               && (canonicalSource.equalsIgnoreCase(DBRefSource.UNIPROT)
 408                       || canonicalSource.equalsIgnoreCase(DBRefSource.PDB)))
 409       {
 410         return dbRef;
 411       }
 412     }
 413     throw new SiftsException("Could not get source DB Ref");
 414   }
 415
 416   /**
 417    * Check that the DBRef Entry is properly populated and is available in this
 418    * SiftClient instance
 419    *
 420    * @param entry
 421    *          - DBRefEntry to validate
 422    * @return true validation is successful otherwise false is returned.
 423    */
 424   boolean isValidDBRefEntry(DBRefEntryI entry)
 425   {
 426     return entry != null && entry.getAccessionId() != null
 427             && isFoundInSiftsEntry(entry.getAccessionId());
 428   }
 429
 430   @Override
 431   public HashSet<String> getAllMappingAccession()
 432   {
 433     HashSet<String> accessions = new HashSet<String>();
 434     List<Entity> entities = siftsEntry.getEntity();
 435     for (Entity entity : entities)
 436     {
 437       List<Segment> segments = entity.getSegment();
 438       for (Segment segment : segments)
 439       {
 440         List<MapRegion> mapRegions = segment.getListMapRegion()
 441                 .getMapRegion();
 442         for (MapRegion mapRegion : mapRegions)
 443         {
 444           accessions.add(mapRegion.getDb().getDbAccessionId()
 445                   .toLowerCase(Locale.ROOT));
 446         }
 447       }
 448     }
 449     return accessions;
 450   }
 451
 452   @Override
 453   public StructureMapping getSiftsStructureMapping(SequenceI seq,
 454           String pdbFile, String chain) throws SiftsException
 455   {
 456     SequenceI aseq = seq;
 457     while (seq.getDatasetSequence() != null)
 458     {
 459       seq = seq.getDatasetSequence();
 460     }
 461     structId = (chain == null) ? pdbId : pdbId + "|" + chain;
 462     System.out.println("Getting SIFTS mapping for " + structId + ": seq "
 463             + seq.getName());
 464
 465     final StringBuilder mappingDetails = new StringBuilder(128);
 466     PrintStream ps = new PrintStream(System.out)
 467     {
 468       @Override
 469       public void print(String x)
 470       {
 471         mappingDetails.append(x);
 472       }
 473
 474       @Override
 475       public void println()
 476       {
 477         mappingDetails.append(NEWLINE);
 478       }
 479     };
 480     HashMap<Integer, int[]> mapping = getGreedyMapping(chain, seq, ps);
 481
 482     String mappingOutput = mappingDetails.toString();
 483     StructureMapping siftsMapping = new StructureMapping(aseq, pdbFile,
 484             pdbId, chain, mapping, mappingOutput, seqFromPdbMapping);
 485
 486     return siftsMapping;
 487   }
 488
 489   @Override
 490   public HashMap<Integer, int[]> getGreedyMapping(String entityId,
 491           SequenceI seq, java.io.PrintStream os) throws SiftsException
 492   {
 493     List<Integer> omitNonObserved = new ArrayList<>();
 494     int nonObservedShiftIndex = 0, pdbeNonObserved = 0;
 495     // System.out.println("Generating mappings for : " + entityId);
 496     Entity entity = null;
 497     entity = getEntityById(entityId);
 498     String originalSeq = AlignSeq.extractGaps(
 499             jalview.util.Comparison.GapChars, seq.getSequenceAsString());
 500     HashMap<Integer, int[]> mapping = new HashMap<Integer, int[]>();
 501     DBRefEntryI sourceDBRef;
 502     sourceDBRef = getValidSourceDBRef(seq);
 503     // TODO ensure sequence start/end is in the same coordinate system and
 504     // consistent with the choosen sourceDBRef
 505
 506     // set sequence coordinate system - default value is UniProt
 507     if (sourceDBRef.getSource().equalsIgnoreCase(DBRefSource.PDB))
 508     {
 509       seqCoordSys = CoordinateSys.PDB;
 510     }
 511
 512     HashSet<String> dbRefAccessionIdsString = new HashSet<String>();
 513     for (DBRefEntry dbref : seq.getDBRefs())
 514     {
 515       dbRefAccessionIdsString
 516               .add(dbref.getAccessionId().toLowerCase(Locale.ROOT));
 517     }
 518     dbRefAccessionIdsString
 519             .add(sourceDBRef.getAccessionId().toLowerCase(Locale.ROOT));
 520
 521     curDBRefAccessionIdsString = dbRefAccessionIdsString;
 522     curSourceDBRef = sourceDBRef.getAccessionId();
 523
 524     TreeMap<Integer, String> resNumMap = new TreeMap<Integer, String>();
 525     List<Segment> segments = entity.getSegment();
 526     SegmentHelperPojo shp = new SegmentHelperPojo(seq, mapping, resNumMap,
 527             omitNonObserved, nonObservedShiftIndex, pdbeNonObserved);
 528     processSegments(segments, shp);
 529     try
 530     {
 531       populateAtomPositions(entityId, mapping);
 532     } catch (Exception e)
 533     {
 534       e.printStackTrace();
 535     }
 536     if (seqCoordSys == CoordinateSys.UNIPROT)
 537     {
 538       padWithGaps(resNumMap, omitNonObserved);
 539     }
 540     int seqStart = UNASSIGNED;
 541     int seqEnd = UNASSIGNED;
 542     int pdbStart = UNASSIGNED;
 543     int pdbEnd = UNASSIGNED;
 544
 545     if (mapping.isEmpty())
 546     {
 547       throw new SiftsException("SIFTS mapping failed for "+entityId+" and "+seq.getName());
 548     }
 549     // also construct a mapping object between the seq-coord sys and the PDB
 550     // seq's coord sys
 551
 552     Integer[] keys = mapping.keySet().toArray(new Integer[0]);
 553     Arrays.sort(keys);
 554     seqStart = keys[0];
 555     seqEnd = keys[keys.length - 1];
 556     List<int[]> from = new ArrayList<>(), to = new ArrayList<>();
 557     int[] _cfrom = null, _cto = null;
 558     String matchedSeq = originalSeq;
 559     if (seqStart != UNASSIGNED) // fixme! seqStart can map to -1 for a pdb
 560                                 // sequence that starts <-1
 561     {
 562       for (int seqps : keys)
 563       {
 564         int pdbpos = mapping.get(seqps)[PDBE_POS];
 565         if (pdbpos == UNASSIGNED)
 566         {
 567           // not correct - pdbpos might be -1, but leave it for now
 568           continue;
 569         }
 570         if (_cfrom == null || seqps != _cfrom[1] + 1)
 571         {
 572           _cfrom = new int[] { seqps, seqps };
 573           from.add(_cfrom);
 574           _cto = null; // discontinuity
 575         }
 576         else
 577         {
 578           _cfrom[1] = seqps;
 579         }
 580         if (_cto == null || pdbpos != 1 + _cto[1])
 581         {
 582           _cto = new int[] { pdbpos, pdbpos };
 583           to.add(_cto);
 584         }
 585         else
 586         {
 587           _cto[1] = pdbpos;
 588         }
 589       }
 590       _cfrom = new int[from.size() * 2];
 591       _cto = new int[to.size() * 2];
 592       int p = 0;
 593       for (int[] range : from)
 594       {
 595         _cfrom[p++] = range[0];
 596         _cfrom[p++] = range[1];
 597       }
 598       ;
 599       p = 0;
 600       for (int[] range : to)
 601       {
 602         _cto[p++] = range[0];
 603         _cto[p++] = range[1];
 604       }
 605       ;
 606
 607       seqFromPdbMapping = new jalview.datamodel.Mapping(null, _cto, _cfrom,
 608               1, 1);
 609       pdbStart = mapping.get(seqStart)[PDB_RES_POS];
 610       pdbEnd = mapping.get(seqEnd)[PDB_RES_POS];
 611       int orignalSeqStart = seq.getStart();
 612       if (orignalSeqStart >= 1)
 613       {
 614         int subSeqStart = (seqStart >= orignalSeqStart)
 615                 ? seqStart - orignalSeqStart
 616                 : 0;
 617         int subSeqEnd = seqEnd - (orignalSeqStart - 1);
 618         subSeqEnd = originalSeq.length() < subSeqEnd ? originalSeq.length()
 619                 : subSeqEnd;
 620         matchedSeq = originalSeq.substring(subSeqStart, subSeqEnd);
 621       }
 622       else
 623       {
 624         matchedSeq = originalSeq.substring(1, originalSeq.length());
 625       }
 626     }
 627
 628     StringBuilder targetStrucSeqs = new StringBuilder();
 629     for (String res : resNumMap.values())
 630     {
 631       targetStrucSeqs.append(res);
 632     }
 633
 634     if (os != null)
 635     {
 636       MappingOutputPojo mop = new MappingOutputPojo();
 637       mop.setSeqStart(seqStart);
 638       mop.setSeqEnd(seqEnd);
 639       mop.setSeqName(seq.getName());
 640       mop.setSeqResidue(matchedSeq);
 641
 642       mop.setStrStart(pdbStart);
 643       mop.setStrEnd(pdbEnd);
 644       mop.setStrName(structId);
 645       mop.setStrResidue(targetStrucSeqs.toString());
 646
 647       mop.setType("pep");
 648       os.print(getMappingOutput(mop).toString());
 649       os.println();
 650     }
 651     return mapping;
 652   }
 653
 654   void processSegments(List<Segment> segments, SegmentHelperPojo shp)
 655   {
 656     SequenceI seq = shp.getSeq();
 657     HashMap<Integer, int[]> mapping = shp.getMapping();
 658     TreeMap<Integer, String> resNumMap = shp.getResNumMap();
 659     List<Integer> omitNonObserved = shp.getOmitNonObserved();
 660     int nonObservedShiftIndex = shp.getNonObservedShiftIndex();
 661     int pdbeNonObservedCount = shp.getPdbeNonObserved();
 662     int firstPDBResNum = UNASSIGNED;
 663     for (Segment segment : segments)
 664     {
 665       // System.out.println("Mapping segments : " + segment.getSegId() + "\\"s
 666       // + segStartEnd);
 667       List<Residue> residues = segment.getListResidue().getResidue();
 668       for (Residue residue : residues)
 669       {
 670         boolean isObserved = isResidueObserved(residue);
 671         int pdbeIndex = getLeadingIntegerValue(residue.getDbResNum(),
 672                 UNASSIGNED);
 673         int currSeqIndex = UNASSIGNED;
 674         List<CrossRefDb> cRefDbs = residue.getCrossRefDb();
 675         CrossRefDb pdbRefDb = null;
 676         for (CrossRefDb cRefDb : cRefDbs)
 677         {
 678           if (cRefDb.getDbSource().equalsIgnoreCase(DBRefSource.PDB))
 679           {
 680             pdbRefDb = cRefDb;
 681             if (firstPDBResNum == UNASSIGNED)
 682             {
 683               firstPDBResNum = getLeadingIntegerValue(cRefDb.getDbResNum(),
 684                       UNASSIGNED);
 685             }
 686             else
 687             {
 688               if (isObserved)
 689               {
 690                 // after we find the first observed residue we just increment
 691                 firstPDBResNum++;
 692               }
 693             }
 694           }
 695           if (cRefDb.getDbCoordSys().equalsIgnoreCase(seqCoordSys.getName())
 696                   && isAccessionMatched(cRefDb.getDbAccessionId()))
 697           {
 698             currSeqIndex = getLeadingIntegerValue(cRefDb.getDbResNum(),
 699                     UNASSIGNED);
 700             if (pdbRefDb != null)
 701             {
 702               break;// exit loop if pdb and uniprot are already found
 703             }
 704           }
 705         }
 706         if (!isObserved)
 707         {
 708           ++pdbeNonObservedCount;
 709         }
 710         if (seqCoordSys == seqCoordSys.PDB) // FIXME: is seqCoordSys ever PDBe
 711                                             // ???
 712         {
 713           // if the sequence has a primary reference to the PDB, then we are
 714           // dealing with a sequence extracted directly from the PDB. In that
 715           // case, numbering is PDBe - non-observed residues
 716           currSeqIndex = seq.getStart() - 1 + pdbeIndex;
 717         }
 718         if (!isObserved)
 719         {
 720           if (seqCoordSys != CoordinateSys.UNIPROT) // FIXME: PDB or PDBe only
 721                                                     // here
 722           {
 723             // mapping to PDB or PDBe so we need to bookkeep for the
 724             // non-observed
 725             // SEQRES positions
 726             omitNonObserved.add(currSeqIndex);
 727             ++nonObservedShiftIndex;
 728           }
 729         }
 730         if (currSeqIndex == UNASSIGNED)
 731         {
 732           // change in logic - unobserved residues with no currSeqIndex
 733           // corresponding are still counted in both nonObservedShiftIndex and
 734           // pdbeIndex...
 735           continue;
 736         }
 737         // if (currSeqIndex >= seq.getStart() && currSeqIndex <= seqlength) //
 738         // true
 739         // numbering
 740         // is
 741         // not
 742         // up
 743         // to
 744         // seq.getEnd()
 745         {
 746
 747           int resNum = (pdbRefDb == null)
 748                   ? getLeadingIntegerValue(residue.getDbResNum(),
 749                           UNASSIGNED)
 750                   : getLeadingIntegerValue(pdbRefDb.getDbResNum(),
 751                           UNASSIGNED);
 752
 753           if (isObserved)
 754           {
 755             char resCharCode = ResidueProperties
 756                     .getSingleCharacterCode(ResidueProperties
 757                             .getCanonicalAminoAcid(residue.getDbResName()));
 758             resNumMap.put(currSeqIndex, String.valueOf(resCharCode));
 759
 760             int[] mappingcols = new int[] { Integer.valueOf(resNum),
 761                 UNASSIGNED, isObserved ? firstPDBResNum : UNASSIGNED };
 762
 763             mapping.put(currSeqIndex - nonObservedShiftIndex, mappingcols);
 764           }
 765         }
 766       }
 767     }
 768   }
 769
 770   /**
 771    * Get the leading integer part of a string that begins with an integer.
 772    *
 773    * @param input
 774    *          - the string input to process
 775    * @param failValue
 776    *          - value returned if unsuccessful
 777    * @return
 778    */
 779   static int getLeadingIntegerValue(String input, int failValue)
 780   {
 781     if (input == null)
 782     {
 783       return failValue;
 784     }
 785     String[] parts = input.split("(?=\\D)(?<=\\d)");
 786     if (parts != null && parts.length > 0 && parts[0].matches("[0-9]+"))
 787     {
 788       return Integer.valueOf(parts[0]);
 789     }
 790     return failValue;
 791   }
 792
 793   /**
 794    *
 795    * @param chainId
 796    *          Target chain to populate mapping of its atom positions.
 797    * @param mapping
 798    *          Two dimension array of residue index versus atom position
 799    * @throws IllegalArgumentException
 800    *           Thrown if chainId or mapping is null
 801    * @throws SiftsException
 802    */
 803   void populateAtomPositions(String chainId, Map<Integer, int[]> mapping)
 804           throws IllegalArgumentException, SiftsException
 805   {
 806     try
 807     {
 808       PDBChain chain = pdb.findChain(chainId);
 809
 810       if (chain == null || mapping == null)
 811       {
 812         throw new IllegalArgumentException(
 813                 "Chain id or mapping must not be null.");
 814       }
 815       for (int[] map : mapping.values())
 816       {
 817         if (map[PDB_RES_POS] != UNASSIGNED)
 818         {
 819           map[PDB_ATOM_POS] = getAtomIndex(map[PDB_RES_POS], chain.atoms);
 820         }
 821       }
 822     } catch (NullPointerException e)
 823     {
 824       throw new SiftsException(e.getMessage());
 825     } catch (Exception e)
 826     {
 827       throw new SiftsException(e.getMessage());
 828     }
 829   }
 830
 831   /**
 832    *
 833    * @param residueIndex
 834    *          The residue index used for the search
 835    * @param atoms
 836    *          A collection of Atom to search
 837    * @return atom position for the given residue index
 838    */
 839   int getAtomIndex(int residueIndex, Collection<Atom> atoms)
 840   {
 841     if (atoms == null)
 842     {
 843       throw new IllegalArgumentException(
 844               "atoms collection must not be null!");
 845     }
 846     for (Atom atom : atoms)
 847     {
 848       if (atom.resNumber == residueIndex)
 849       {
 850         return atom.atomIndex;
 851       }
 852     }
 853     return UNASSIGNED;
 854   }
 855
 856   /**
 857    * Checks if the residue instance is marked 'Not_observed' or not
 858    *
 859    * @param residue
 860    * @return
 861    */
 862   private boolean isResidueObserved(Residue residue)
 863   {
 864     Set<String> annotations = getResidueAnnotaitons(residue,
 865             ResidueDetailType.ANNOTATION);
 866     if (annotations == null || annotations.isEmpty())
 867     {
 868       return true;
 869     }
 870     for (String annotation : annotations)
 871     {
 872       if (annotation.equalsIgnoreCase(NOT_OBSERVED))
 873       {
 874         return false;
 875       }
 876     }
 877     return true;
 878   }
 879
 880   /**
 881    * Get annotation String for a given residue and annotation type
 882    *
 883    * @param residue
 884    * @param type
 885    * @return
 886    */
 887   private Set<String> getResidueAnnotaitons(Residue residue,
 888           ResidueDetailType type)
 889   {
 890     HashSet<String> foundAnnotations = new HashSet<String>();
 891     List<ResidueDetail> resDetails = residue.getResidueDetail();
 892     for (ResidueDetail resDetail : resDetails)
 893     {
 894       if (resDetail.getProperty().equalsIgnoreCase(type.getCode()))
 895       {
 896         foundAnnotations.add(resDetail.getContent());
 897       }
 898     }
 899     return foundAnnotations;
 900   }
 901
 902   @Override
 903   public boolean isAccessionMatched(String accession)
 904   {
 905     boolean isStrictMatch = true;
 906     return isStrictMatch ? curSourceDBRef.equalsIgnoreCase(accession)
 907             : curDBRefAccessionIdsString
 908                     .contains(accession.toLowerCase(Locale.ROOT));
 909   }
 910
 911   private boolean isFoundInSiftsEntry(String accessionId)
 912   {
 913     Set<String> siftsDBRefs = getAllMappingAccession();
 914     return accessionId != null
 915             && siftsDBRefs.contains(accessionId.toLowerCase(Locale.ROOT));
 916   }
 917
 918   /**
 919    * Pad omitted residue positions in PDB sequence with gaps
 920    *
 921    * @param resNumMap
 922    */
 923   void padWithGaps(Map<Integer, String> resNumMap,
 924           List<Integer> omitNonObserved)
 925   {
 926     if (resNumMap == null || resNumMap.isEmpty())
 927     {
 928       return;
 929     }
 930     Integer[] keys = resNumMap.keySet().toArray(new Integer[0]);
 931     // Arrays.sort(keys);
 932     int firstIndex = keys[0];
 933     int lastIndex = keys[keys.length - 1];
 934     // System.out.println("Min value " + firstIndex);
 935     // System.out.println("Max value " + lastIndex);
 936     for (int x = firstIndex; x <= lastIndex; x++)
 937     {
 938       if (!resNumMap.containsKey(x) && !omitNonObserved.contains(x))
 939       {
 940         resNumMap.put(x, "-");
 941       }
 942     }
 943   }
 944
 945   @Override
 946   public Entity getEntityById(String id) throws SiftsException
 947   {
 948     // Determines an entity to process by performing a heuristic matching of all
 949     // Entities with the given chainId and choosing the best matching Entity
 950     Entity entity = getEntityByMostOptimalMatchedId(id);
 951     if (entity != null)
 952     {
 953       return entity;
 954     }
 955     throw new SiftsException("Entity " + id + " not found");
 956   }
 957
 958   /**
 959    * This method was added because EntityId is NOT always equal to ChainId.
 960    * Hence, it provides the logic to greedily detect the "true" Entity for a
 961    * given chainId where discrepancies exist.
 962    *
 963    * @param chainId
 964    * @return
 965    */
 966   public Entity getEntityByMostOptimalMatchedId(String chainId)
 967   {
 968     // System.out.println("---> advanced greedy entityId matching block
 969     // entered..");
 970     List<Entity> entities = siftsEntry.getEntity();
 971     SiftsEntitySortPojo[] sPojo = new SiftsEntitySortPojo[entities.size()];
 972     int count = 0;
 973     for (Entity entity : entities)
 974     {
 975       sPojo[count] = new SiftsEntitySortPojo();
 976       sPojo[count].entityId = entity.getEntityId();
 977
 978       List<Segment> segments = entity.getSegment();
 979       for (Segment segment : segments)
 980       {
 981         List<Residue> residues = segment.getListResidue().getResidue();
 982         for (Residue residue : residues)
 983         {
 984           List<CrossRefDb> cRefDbs = residue.getCrossRefDb();
 985           for (CrossRefDb cRefDb : cRefDbs)
 986           {
 987             if (!cRefDb.getDbSource().equalsIgnoreCase("PDB"))
 988             {
 989               continue;
 990             }
 991             ++sPojo[count].resCount;
 992             if (cRefDb.getDbChainId().equalsIgnoreCase(chainId))
 993             {
 994               ++sPojo[count].chainIdFreq;
 995             }
 996           }
 997         }
 998       }
 999       sPojo[count].pid = (100 * sPojo[count].chainIdFreq)
1000               / sPojo[count].resCount;
1001       ++count;
1002     }
1003     Arrays.sort(sPojo, Collections.reverseOrder());
1004     // System.out.println("highest matched entity : " + sPojo[0].entityId);
1005     // System.out.println("highest matched pid : " + sPojo[0].pid);
1006
1007     if (sPojo[0].entityId != null)
1008     {
1009       if (sPojo[0].pid < 1)
1010       {
1011         return null;
1012       }
1013       for (Entity entity : entities)
1014       {
1015         if (!entity.getEntityId().equalsIgnoreCase(sPojo[0].entityId))
1016         {
1017           continue;
1018         }
1019         return entity;
1020       }
1021     }
1022     return null;
1023   }
1024
1025   private class SiftsEntitySortPojo
1026           implements Comparable<SiftsEntitySortPojo>
1027   {
1028     public String entityId;
1029
1030     public int chainIdFreq;
1031
1032     public int pid;
1033
1034     public int resCount;
1035
1036     @Override
1037     public int compareTo(SiftsEntitySortPojo o)
1038     {
1039       return this.pid - o.pid;
1040     }
1041   }
1042
1043   private class SegmentHelperPojo
1044   {
1045     private SequenceI seq;
1046
1047     private HashMap<Integer, int[]> mapping;
1048
1049     private TreeMap<Integer, String> resNumMap;
1050
1051     private List<Integer> omitNonObserved;
1052
1053     private int nonObservedShiftIndex;
1054
1055     /**
1056      * count of number of 'not observed' positions in the PDB record's SEQRES
1057      * (total number of residues with coordinates == length(SEQRES) -
1058      * pdbeNonObserved
1059      */
1060     private int pdbeNonObserved;
1061
1062     public SegmentHelperPojo(SequenceI seq, HashMap<Integer, int[]> mapping,
1063             TreeMap<Integer, String> resNumMap,
1064             List<Integer> omitNonObserved, int nonObservedShiftIndex,
1065             int pdbeNonObserved)
1066     {
1067       setSeq(seq);
1068       setMapping(mapping);
1069       setResNumMap(resNumMap);
1070       setOmitNonObserved(omitNonObserved);
1071       setNonObservedShiftIndex(nonObservedShiftIndex);
1072       setPdbeNonObserved(pdbeNonObserved);
1073
1074     }
1075
1076     public void setPdbeNonObserved(int pdbeNonObserved2)
1077     {
1078       this.pdbeNonObserved = pdbeNonObserved2;
1079     }
1080
1081     public int getPdbeNonObserved()
1082     {
1083       return pdbeNonObserved;
1084     }
1085
1086     public SequenceI getSeq()
1087     {
1088       return seq;
1089     }
1090
1091     public void setSeq(SequenceI seq)
1092     {
1093       this.seq = seq;
1094     }
1095
1096     public HashMap<Integer, int[]> getMapping()
1097     {
1098       return mapping;
1099     }
1100
1101     public void setMapping(HashMap<Integer, int[]> mapping)
1102     {
1103       this.mapping = mapping;
1104     }
1105
1106     public TreeMap<Integer, String> getResNumMap()
1107     {
1108       return resNumMap;
1109     }
1110
1111     public void setResNumMap(TreeMap<Integer, String> resNumMap)
1112     {
1113       this.resNumMap = resNumMap;
1114     }
1115
1116     public List<Integer> getOmitNonObserved()
1117     {
1118       return omitNonObserved;
1119     }
1120
1121     public void setOmitNonObserved(List<Integer> omitNonObserved)
1122     {
1123       this.omitNonObserved = omitNonObserved;
1124     }
1125
1126     public int getNonObservedShiftIndex()
1127     {
1128       return nonObservedShiftIndex;
1129     }
1130
1131     public void setNonObservedShiftIndex(int nonObservedShiftIndex)
1132     {
1133       this.nonObservedShiftIndex = nonObservedShiftIndex;
1134     }
1135
1136   }
1137
1138   @Override
1139   public StringBuilder getMappingOutput(MappingOutputPojo mp)
1140           throws SiftsException
1141   {
1142     String seqRes = mp.getSeqResidue();
1143     String seqName = mp.getSeqName();
1144     int sStart = mp.getSeqStart();
1145     int sEnd = mp.getSeqEnd();
1146
1147     String strRes = mp.getStrResidue();
1148     String strName = mp.getStrName();
1149     int pdbStart = mp.getStrStart();
1150     int pdbEnd = mp.getStrEnd();
1151
1152     String type = mp.getType();
1153
1154     int maxid = (seqName.length() >= strName.length()) ? seqName.length()
1155             : strName.length();
1156     int len = 72 - maxid - 1;
1157
1158     int nochunks = ((seqRes.length()) / len)
1159             + ((seqRes.length()) % len > 0 ? 1 : 0);
1160     // output mappings
1161     StringBuilder output = new StringBuilder(512);
1162     output.append(NEWLINE);
1163     output.append("Sequence \u27f7 Structure mapping details")
1164             .append(NEWLINE);
1165     output.append("Method: SIFTS");
1166     output.append(NEWLINE).append(NEWLINE);
1167
1168     output.append(new Format("%" + maxid + "s").form(seqName));
1169     output.append(" :  ");
1170     output.append(String.valueOf(sStart));
1171     output.append(" - ");
1172     output.append(String.valueOf(sEnd));
1173     output.append(" Maps to ");
1174     output.append(NEWLINE);
1175     output.append(new Format("%" + maxid + "s").form(structId));
1176     output.append(" :  ");
1177     output.append(String.valueOf(pdbStart));
1178     output.append(" - ");
1179     output.append(String.valueOf(pdbEnd));
1180     output.append(NEWLINE).append(NEWLINE);
1181
1182     ScoreMatrix pam250 = ScoreModels.getInstance().getPam250();
1183     int matchedSeqCount = 0;
1184     for (int j = 0; j < nochunks; j++)
1185     {
1186       // Print the first aligned sequence
1187       output.append(new Format("%" + (maxid) + "s").form(seqName))
1188               .append(" ");
1189
1190       for (int i = 0; i < len; i++)
1191       {
1192         if ((i + (j * len)) < seqRes.length())
1193         {
1194           output.append(seqRes.charAt(i + (j * len)));
1195         }
1196       }
1197
1198       output.append(NEWLINE);
1199       output.append(new Format("%" + (maxid) + "s").form(" ")).append(" ");
1200
1201       /*
1202        * Print out the match symbols:
1203        * | for exact match (ignoring case)
1204        * . if PAM250 score is positive
1205        * else a space
1206        */
1207       for (int i = 0; i < len; i++)
1208       {
1209         try
1210         {
1211           if ((i + (j * len)) < seqRes.length())
1212           {
1213             char c1 = seqRes.charAt(i + (j * len));
1214             char c2 = strRes.charAt(i + (j * len));
1215             boolean sameChar = Comparison.isSameResidue(c1, c2, false);
1216             if (sameChar && !Comparison.isGap(c1))
1217             {
1218               matchedSeqCount++;
1219               output.append("|");
1220             }
1221             else if (type.equals("pep"))
1222             {
1223               if (pam250.getPairwiseScore(c1, c2) > 0)
1224               {
1225                 output.append(".");
1226               }
1227               else
1228               {
1229                 output.append(" ");
1230               }
1231             }
1232             else
1233             {
1234               output.append(" ");
1235             }
1236           }
1237         } catch (IndexOutOfBoundsException e)
1238         {
1239           continue;
1240         }
1241       }
1242       // Now print the second aligned sequence
1243       output = output.append(NEWLINE);
1244       output = output.append(new Format("%" + (maxid) + "s").form(strName))
1245               .append(" ");
1246       for (int i = 0; i < len; i++)
1247       {
1248         if ((i + (j * len)) < strRes.length())
1249         {
1250           output.append(strRes.charAt(i + (j * len)));
1251         }
1252       }
1253       output.append(NEWLINE).append(NEWLINE);
1254     }
1255     float pid = (float) matchedSeqCount / seqRes.length() * 100;
1256     if (pid < SiftsSettings.getFailSafePIDThreshold())
1257     {
1258       throw new SiftsException(">>> Low PID detected for SIFTs mapping...");
1259     }
1260     output.append("Length of alignment = " + seqRes.length())
1261             .append(NEWLINE);
1262     output.append(new Format("Percentage ID = %2.2f").form(pid));
1263     return output;
1264   }
1265
1266   @Override
1267   public int getEntityCount()
1268   {
1269     return siftsEntry.getEntity().size();
1270   }
1271
1272   @Override
1273   public String getDbAccessionId()
1274   {
1275     return siftsEntry.getDbAccessionId();
1276   }
1277
1278   @Override
1279   public String getDbCoordSys()
1280   {
1281     return siftsEntry.getDbCoordSys();
1282   }
1283
1284   @Override
1285   public String getDbSource()
1286   {
1287     return siftsEntry.getDbSource();
1288   }
1289
1290   @Override
1291   public String getDbVersion()
1292   {
1293     return siftsEntry.getDbVersion();
1294   }
1295
1296   public static void setMockSiftsFile(File file)
1297   {
1298     mockSiftsFile = file;
1299   }
1300
1301 }