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