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