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