c64dac1fd2cb75249fbc98958227beb1c50ad296
[jalview.git] / src / jalview / ext / jmol / JmolParser.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.ext.jmol;
22
23 import java.io.File;
24 import java.io.IOException;
25 import java.util.ArrayList;
26 import java.util.HashMap;
27 import java.util.List;
28 import java.util.Locale;
29 import java.util.Map;
30 import java.util.Vector;
31
32 import org.jmol.api.JmolStatusListener;
33 import org.jmol.api.JmolViewer;
34 import org.jmol.c.CBK;
35 import org.jmol.c.STR;
36 import org.jmol.modelset.ModelSet;
37 import org.jmol.viewer.Viewer;
38
39 import com.stevesoft.pat.Regex;
40
41 import jalview.bin.Console;
42 import jalview.datamodel.Alignment;
43 import jalview.datamodel.AlignmentAnnotation;
44 import jalview.datamodel.Annotation;
45 import jalview.datamodel.PDBEntry;
46 import jalview.datamodel.SequenceI;
47 import jalview.datamodel.annotations.AlphaFoldAnnotationRowBuilder;
48 import jalview.datamodel.annotations.AnnotationRowBuilder;
49 import jalview.io.DataSourceType;
50 import jalview.io.FileParse;
51 import jalview.io.StructureFile;
52 import jalview.schemes.ResidueProperties;
53 import jalview.structure.StructureImportSettings;
54 import jalview.util.Format;
55 import jalview.util.MessageManager;
56 import jalview.ws.dbsources.EBIAlfaFold;
57 import mc_view.Atom;
58 import mc_view.PDBChain;
59 import mc_view.Residue;
60
61 /**
62  * Import and process files with Jmol for file like PDB, mmCIF
63  * 
64  * @author jprocter
65  * 
66  */
67 public class JmolParser extends StructureFile implements JmolStatusListener
68 {
69   Viewer viewer = null;
70
71   public JmolParser(boolean immediate, Object inFile,
72           DataSourceType sourceType) throws IOException
73   {
74     // BH 2018 File or String for filename
75     super(immediate, inFile, sourceType);
76
77   }
78
79   public JmolParser(Object inFile, DataSourceType sourceType)
80           throws IOException
81   {
82     this(inFile, sourceType, null);
83   }
84
85   public JmolParser(Object inFile, DataSourceType sourceType,
86           StructureImportSettings.TFType tempfacType) throws IOException
87   {
88     super(inFile, sourceType, tempfacType);
89   }
90
91   public JmolParser(FileParse fp, boolean doXferSettings) throws IOException
92   {
93     super(fp, doXferSettings);
94   }
95   public JmolParser(FileParse fp) throws IOException
96   {
97     super(fp);
98   }
99
100   public JmolParser()
101   {
102   }
103
104   /**
105    * Calls the Jmol library to parse the PDB/mmCIF file, and then inspects the
106    * resulting object model to generate Jalview-style sequences, with secondary
107    * structure annotation added where available (i.e. where it has been computed
108    * by Jmol using DSSP).
109    * 
110    * @see jalview.io.AlignFile#parse()
111    */
112   @Override
113   public void parse() throws IOException
114   {
115     parse(true);
116   }
117
118   @Override
119   public void parse(boolean doXferSettings) throws IOException
120   {
121     setChains(new Vector<PDBChain>());
122     Viewer jmolModel = getJmolData();
123     jmolModel.openReader(getDataName(), getDataName(), getReader());
124     waitForScript(jmolModel);
125
126     /*
127      * Convert one or more Jmol Model objects to Jalview sequences
128      */
129     if (jmolModel.ms.mc > 0)
130     {
131       // ideally we do this
132       // try
133       // {
134       // setStructureFileType(jmolModel.evalString("show _fileType"));
135       // } catch (Exception q)
136       // {
137       // }
138       // ;
139       // instead, we distinguish .cif from non-.cif by filename
140       setStructureFileType(
141               getDataName().toLowerCase(Locale.ROOT).endsWith(".cif")
142                       ? PDBEntry.Type.MMCIF.toString()
143                       : "PDB");
144
145       transformJmolModelToJalview(jmolModel.ms, doXferSettings);
146     }
147   }
148
149   /**
150    * create a headless jmol instance for dataprocessing
151    * 
152    * @return
153    */
154   private Viewer getJmolData()
155   {
156     if (viewer == null)
157     {
158       try
159       {
160         /*
161          * params -o (output to sysout) -n (nodisplay) -x (exit when finished)
162          * see http://wiki.jmol.org/index.php/Jmol_Application
163          */
164
165         viewer = JalviewJmolBinding.getJmolData(this);
166         // ensure the 'new' (DSSP) not 'old' (Ramachandran) SS method is used
167         viewer.setBooleanProperty("defaultStructureDSSP", true);
168       } catch (ClassCastException x)
169       {
170         throw new Error(MessageManager.formatMessage(
171                 "error.jmol_version_not_compatible_with_jalview_version",
172                 new String[]
173                 { JmolViewer.getJmolVersion() }), x);
174       }
175     }
176     return viewer;
177   }
178
179   public static Regex getNewAlphafoldValidator()
180   {
181     Regex validator = new Regex("(AF-[A-Z]+[0-9]+[A-Z0-9]+-F1)");
182     validator.setIgnoreCase(true);
183     return validator;
184   }
185
186   PDBEntry.Type jmolFiletype = null;
187
188   /**
189    * resolve a jmol filetype string and update the jmolFiletype field
190    * accordingly
191    * 
192    * @param jmolIdentifiedFileType
193    * @return true if filetype was identified as MMCIF, PDB
194    */
195   public boolean updateFileType(String jmolIdentifiedFileType)
196   {
197     if (jmolIdentifiedFileType == null
198             || jmolIdentifiedFileType.trim().equals(""))
199     {
200       return false;
201     }
202     if ("mmcif".equalsIgnoreCase(jmolIdentifiedFileType))
203     {
204       jmolFiletype = PDBEntry.Type.MMCIF;
205       return true;
206     }
207     if ("pdb".equalsIgnoreCase(jmolIdentifiedFileType))
208     {
209       jmolFiletype = PDBEntry.Type.PDB;
210       return true;
211     }
212     return false;
213   }
214
215   public void transformJmolModelToJalview(ModelSet ms,
216           boolean localDoXferSettings) throws IOException
217   {
218     try
219     {
220       Regex alphaFold = getNewAlphafoldValidator();
221       String lastID = "";
222       List<SequenceI> rna = new ArrayList<SequenceI>();
223       List<SequenceI> prot = new ArrayList<SequenceI>();
224       PDBChain tmpchain;
225       String pdbId = (String) ms.getInfo(0, "title");
226       boolean isMMCIF = false;
227       String jmolFileType_String = (String) ms.getInfo(0, "fileType");
228       if (updateFileType(jmolFileType_String))
229       {
230         setStructureFileType(jmolFiletype.toString());
231       }
232
233       isMMCIF = PDBEntry.Type.MMCIF.equals(jmolFiletype);
234
235       if (pdbId == null)
236       {
237         setId(safeName(getDataName()));
238         setPDBIdAvailable(false);
239       }
240       else
241       {
242         setId(pdbId);
243         setPDBIdAvailable(true);
244         setAlphafoldModel(alphaFold.search(pdbId) && isMMCIF);
245       }
246       List<Atom> significantAtoms = convertSignificantAtoms(ms);
247       for (Atom tmpatom : significantAtoms)
248       {
249         if (tmpatom.resNumIns.trim().equals(lastID))
250         {
251           // phosphorylated protein - seen both CA and P..
252           continue;
253         }
254         tmpchain = findChain(tmpatom.chain);
255         if (tmpchain != null)
256         {
257           tmpchain.atoms.addElement(tmpatom);
258         }
259         else
260         {
261           AnnotationRowBuilder builder = null;
262           if (isAlphafoldModel()
263                   || getTemperatureFactorType() == StructureImportSettings.TFType.PLDDT)
264           {
265             builder = new AlphaFoldAnnotationRowBuilder();
266           }
267
268           tmpchain = new PDBChain(getId(), tmpatom.chain, builder);
269           getChains().add(tmpchain);
270           tmpchain.atoms.addElement(tmpatom);
271         }
272         lastID = tmpatom.resNumIns.trim();
273       }
274       if (isParseImmediately() && localDoXferSettings)
275       {
276         // configure parsing settings from the static singleton
277         xferSettings();
278       }
279
280       makeResidueList();
281       makeCaBondList();
282
283       for (PDBChain chain : getChains())
284       {
285         SequenceI chainseq = postProcessChain(chain);
286         if (isRNA(chainseq))
287         {
288           rna.add(chainseq);
289         }
290         else
291         {
292           prot.add(chainseq);
293         }
294
295         // look at local setting for adding secondary tructure
296         if (predictSecondaryStructure)
297         {
298           createAnnotation(chainseq, chain, ms.at);
299         }
300       }
301       // if Alphafold, fetch the PAE matrix if doesn't already have one
302       if (isAlphafoldModel() && !hasPAEMatrix())
303       {
304         try
305         {
306           Console.info("Retrieving PAE for " + pdbId);
307           File paeFile = EBIAlfaFold.fetchAlphaFoldPAE(pdbId, null);
308           this.setPAEMatrix(paeFile.getAbsolutePath());
309         } catch (Throwable t)
310         {
311           Console.error("Couldn't get the pAE for " + pdbId, t);
312         }
313       }
314       // add a PAEMatrix if set (either by above or otherwise)
315       if (hasPAEMatrix())
316       {
317         try
318         {
319           Alignment al = new Alignment(prot.toArray(new SequenceI[0]));
320           EBIAlfaFold.addAlphaFoldPAE(al, new File(this.getPAEMatrix()), 0,
321                   null, false, false, null);
322
323           if (al.getAlignmentAnnotation() != null)
324           {
325             for (AlignmentAnnotation alann : al.getAlignmentAnnotation())
326             {
327               annotations.add(alann);
328             }
329           }
330         } catch (Throwable ff)
331         {
332           Console.error("Couldn't import PAE Matrix from " + getPAEMatrix(),
333                   ff);
334           warningMessage += "Couldn't import PAE Matrix"
335                   + getNewlineString() + ff.getLocalizedMessage()
336                   + getNewlineString();
337         }
338       }
339     } catch (OutOfMemoryError er)
340     {
341       System.out.println(
342               "OUT OF MEMORY LOADING TRANSFORMING JMOL MODEL TO JALVIEW MODEL");
343       throw new IOException(MessageManager
344               .getString("exception.outofmemory_loading_mmcif_file"));
345     }
346   }
347
348   private List<Atom> convertSignificantAtoms(ModelSet ms)
349   {
350     List<Atom> significantAtoms = new ArrayList<Atom>();
351     HashMap<String, org.jmol.modelset.Atom> chainTerMap = new HashMap<String, org.jmol.modelset.Atom>();
352     org.jmol.modelset.Atom prevAtom = null;
353     for (org.jmol.modelset.Atom atom : ms.at)
354     {
355       if (atom.getAtomName().equalsIgnoreCase("CA")
356               || atom.getAtomName().equalsIgnoreCase("P"))
357       {
358         if (!atomValidated(atom, prevAtom, chainTerMap))
359         {
360           continue;
361         }
362         Atom curAtom = new Atom(atom.x, atom.y, atom.z);
363         curAtom.atomIndex = atom.getIndex();
364         curAtom.chain = atom.getChainIDStr();
365         curAtom.insCode = atom.group.getInsertionCode() == '\000' ? ' '
366                 : atom.group.getInsertionCode();
367         curAtom.name = atom.getAtomName();
368         curAtom.number = atom.getAtomNumber();
369         curAtom.resName = atom.getGroup3(true);
370         curAtom.resNumber = atom.getResno();
371         curAtom.occupancy = ms.occupancies != null
372                 ? ms.occupancies[atom.getIndex()]
373                 : Float.valueOf(atom.getOccupancy100());
374         String fmt = new Format("%4i").form(curAtom.resNumber);
375         curAtom.resNumIns = (fmt + curAtom.insCode);
376         curAtom.tfactor = atom.getBfactor100() / 100f;
377         curAtom.type = 0;
378         // significantAtoms.add(curAtom);
379         // ignore atoms from subsequent models
380         if (!significantAtoms.contains(curAtom))
381         {
382           significantAtoms.add(curAtom);
383         }
384         prevAtom = atom;
385       }
386     }
387     return significantAtoms;
388   }
389
390   private boolean atomValidated(org.jmol.modelset.Atom curAtom,
391           org.jmol.modelset.Atom prevAtom,
392           HashMap<String, org.jmol.modelset.Atom> chainTerMap)
393   {
394     // System.out.println("Atom: " + curAtom.getAtomNumber()
395     // + " Last atom index " + curAtom.group.lastAtomIndex);
396     if (chainTerMap == null || prevAtom == null)
397     {
398       return true;
399     }
400     String curAtomChId = curAtom.getChainIDStr();
401     String prevAtomChId = prevAtom.getChainIDStr();
402     // new chain encoutered
403     if (!prevAtomChId.equals(curAtomChId))
404     {
405       // On chain switch add previous chain termination to xTerMap if not exists
406       if (!chainTerMap.containsKey(prevAtomChId))
407       {
408         chainTerMap.put(prevAtomChId, prevAtom);
409       }
410       // if current atom belongs to an already terminated chain and the resNum
411       // diff < 5 then mark as valid and update termination Atom
412       if (chainTerMap.containsKey(curAtomChId))
413       {
414         if (curAtom.getResno() < chainTerMap.get(curAtomChId).getResno())
415         {
416           return false;
417         }
418         if ((curAtom.getResno()
419                 - chainTerMap.get(curAtomChId).getResno()) < 5)
420         {
421           chainTerMap.put(curAtomChId, curAtom);
422           return true;
423         }
424         return false;
425       }
426     }
427     // atom with previously terminated chain encountered
428     else if (chainTerMap.containsKey(curAtomChId))
429     {
430       if (curAtom.getResno() < chainTerMap.get(curAtomChId).getResno())
431       {
432         return false;
433       }
434       if ((curAtom.getResno()
435               - chainTerMap.get(curAtomChId).getResno()) < 5)
436       {
437         chainTerMap.put(curAtomChId, curAtom);
438         return true;
439       }
440       return false;
441     }
442     // HETATM with resNum jump > 2
443     return !(curAtom.isHetero()
444             && ((curAtom.getResno() - prevAtom.getResno()) > 2));
445   }
446
447   private void createAnnotation(SequenceI sequence, PDBChain chain,
448           org.jmol.modelset.Atom[] jmolAtoms)
449   {
450     char[] secstr = new char[sequence.getLength()];
451     char[] secstrcode = new char[sequence.getLength()];
452
453     // Ensure Residue size equals Seq size
454     if (chain.residues.size() != sequence.getLength())
455     {
456       return;
457     }
458     int annotIndex = 0;
459     for (Residue residue : chain.residues)
460     {
461       Atom repAtom = residue.getAtoms().get(0);
462       STR proteinStructureSubType = jmolAtoms[repAtom.atomIndex].group
463               .getProteinStructureSubType();
464       setSecondaryStructure(proteinStructureSubType, annotIndex, secstr,
465               secstrcode);
466       ++annotIndex;
467     }
468     addSecondaryStructureAnnotation(chain.pdbid, sequence, secstr,
469             secstrcode, chain.id, sequence.getStart());
470   }
471
472   /**
473    * Helper method that adds an AlignmentAnnotation for secondary structure to
474    * the sequence, provided at least one secondary structure assignment has been
475    * made
476    * 
477    * @param modelTitle
478    * @param seq
479    * @param secstr
480    * @param secstrcode
481    * @param chainId
482    * @param firstResNum
483    * @return
484    */
485   protected void addSecondaryStructureAnnotation(String modelTitle,
486           SequenceI sq, char[] secstr, char[] secstrcode, String chainId,
487           int firstResNum)
488   {
489     int length = sq.getLength();
490     boolean ssFound = false;
491     Annotation asecstr[] = new Annotation[length + firstResNum - 1];
492     for (int p = 0; p < length; p++)
493     {
494       if (secstr[p] >= 'A' && secstr[p] <= 'z')
495       {
496         try
497         {
498           asecstr[p] = new Annotation(null, null, secstrcode[p], Float.NaN);
499           ssFound = true;
500         } catch (Exception e)
501         {
502           // e.printStackTrace();
503         }
504       }
505     }
506
507     if (ssFound)
508     {
509       String mt = modelTitle == null ? getDataName() : modelTitle;
510       mt += chainId;
511       AlignmentAnnotation ann = new AlignmentAnnotation(
512               "Secondary Structure", "Secondary Structure for " + mt,
513               asecstr);
514       ann.belowAlignment = true;
515       ann.visible = true;
516       ann.autoCalculated = false;
517       ann.setCalcId(getClass().getName());
518       ann.adjustForAlignment();
519       ann.validateRangeAndDisplay();
520       annotations.add(ann);
521       sq.addAlignmentAnnotation(ann);
522     }
523   }
524
525   private void waitForScript(Viewer jmd)
526   {
527     while (jmd.isScriptExecuting())
528     {
529       try
530       {
531         Thread.sleep(50);
532
533       } catch (InterruptedException x)
534       {
535       }
536     }
537   }
538
539   /**
540    * Convert Jmol's secondary structure code to Jalview's, and stored it in the
541    * secondary structure arrays at the given sequence position
542    * 
543    * @param proteinStructureSubType
544    * @param pos
545    * @param secstr
546    * @param secstrcode
547    */
548   protected void setSecondaryStructure(STR proteinStructureSubType, int pos,
549           char[] secstr, char[] secstrcode)
550   {
551     switch (proteinStructureSubType)
552     {
553     case HELIX310:
554       secstr[pos] = '3';
555       break;
556     case HELIX:
557     case HELIXALPHA:
558       secstr[pos] = 'H';
559       break;
560     case HELIXPI:
561       secstr[pos] = 'P';
562       break;
563     case SHEET:
564       secstr[pos] = 'E';
565       break;
566     default:
567       secstr[pos] = 0;
568     }
569
570     switch (proteinStructureSubType)
571     {
572     case HELIX310:
573     case HELIXALPHA:
574     case HELIXPI:
575     case HELIX:
576       secstrcode[pos] = 'H';
577       break;
578     case SHEET:
579       secstrcode[pos] = 'E';
580       break;
581     default:
582       secstrcode[pos] = 0;
583     }
584   }
585
586   /**
587    * Convert any non-standard peptide codes to their standard code table
588    * equivalent. (Initial version only does Selenomethionine MSE->MET.)
589    * 
590    * @param threeLetterCode
591    * @param seq
592    * @param pos
593    */
594   protected void replaceNonCanonicalResidue(String threeLetterCode,
595           char[] seq, int pos)
596   {
597     String canonical = ResidueProperties
598             .getCanonicalAminoAcid(threeLetterCode);
599     if (canonical != null && !canonical.equalsIgnoreCase(threeLetterCode))
600     {
601       seq[pos] = ResidueProperties.getSingleCharacterCode(canonical);
602     }
603   }
604
605   /**
606    * Not implemented - returns null
607    */
608   @Override
609   public String print(SequenceI[] seqs, boolean jvSuffix)
610   {
611     return null;
612   }
613
614   /**
615    * Not implemented
616    */
617   @Override
618   public void setCallbackFunction(String callbackType,
619           String callbackFunction)
620   {
621   }
622
623   @Override
624   public void notifyCallback(CBK cbType, Object[] data)
625   {
626     String strInfo = (data == null || data[1] == null ? null
627             : data[1].toString());
628     switch (cbType)
629     {
630     case ECHO:
631       sendConsoleEcho(strInfo);
632       break;
633     case SCRIPT:
634       notifyScriptTermination((String) data[2],
635               ((Integer) data[3]).intValue());
636       break;
637     case MEASURE:
638       String mystatus = (String) data[3];
639       if (mystatus.indexOf("Picked") >= 0
640               || mystatus.indexOf("Sequence") >= 0)
641       {
642         // Picking mode
643         sendConsoleMessage(strInfo);
644       }
645       else if (mystatus.indexOf("Completed") >= 0)
646       {
647         sendConsoleEcho(strInfo.substring(strInfo.lastIndexOf(",") + 2,
648                 strInfo.length() - 1));
649       }
650       break;
651     case MESSAGE:
652       sendConsoleMessage(data == null ? null : strInfo);
653       break;
654     case PICK:
655       sendConsoleMessage(strInfo);
656       break;
657     default:
658       break;
659     }
660   }
661
662   String lastConsoleEcho = "";
663
664   private void sendConsoleEcho(String string)
665   {
666     lastConsoleEcho += string;
667     lastConsoleEcho += "\n";
668   }
669
670   String lastConsoleMessage = "";
671
672   private void sendConsoleMessage(String string)
673   {
674     lastConsoleMessage += string;
675     lastConsoleMessage += "\n";
676   }
677
678   int lastScriptTermination = -1;
679
680   String lastScriptMessage = "";
681
682   private void notifyScriptTermination(String string, int intValue)
683   {
684     lastScriptMessage += string;
685     lastScriptMessage += "\n";
686     lastScriptTermination = intValue;
687   }
688
689   @Override
690   public boolean notifyEnabled(CBK callbackPick)
691   {
692     switch (callbackPick)
693     {
694     case MESSAGE:
695     case SCRIPT:
696     case ECHO:
697     case LOADSTRUCT:
698     case ERROR:
699       return true;
700     default:
701       return false;
702     }
703   }
704
705   /**
706    * Not implemented - returns null
707    */
708   @Override
709   public String eval(String strEval)
710   {
711     return null;
712   }
713
714   /**
715    * Not implemented - returns null
716    */
717   @Override
718   public float[][] functionXY(String functionName, int x, int y)
719   {
720     return null;
721   }
722
723   /**
724    * Not implemented - returns null
725    */
726   @Override
727   public float[][][] functionXYZ(String functionName, int nx, int ny,
728           int nz)
729   {
730     return null;
731   }
732
733   /**
734    * Not implemented - returns null
735    */
736   @Override
737   public String createImage(String fileName, String imageType,
738           Object text_or_bytes, int quality)
739   {
740     return null;
741   }
742
743   /**
744    * Not implemented - returns null
745    */
746   @Override
747   public Map<String, Object> getRegistryInfo()
748   {
749     return null;
750   }
751
752   /**
753    * Not implemented
754    */
755   @Override
756   public void showUrl(String url)
757   {
758   }
759
760   /**
761    * Not implemented - returns null
762    */
763   @Override
764   public int[] resizeInnerPanel(String data)
765   {
766     return null;
767   }
768
769   @Override
770   public Map<String, Object> getJSpecViewProperty(String arg0)
771   {
772     return null;
773   }
774
775   public boolean isPredictSecondaryStructure()
776   {
777     return predictSecondaryStructure;
778   }
779
780   public void setPredictSecondaryStructure(
781           boolean predictSecondaryStructure)
782   {
783     this.predictSecondaryStructure = predictSecondaryStructure;
784   }
785
786   public boolean isVisibleChainAnnotation()
787   {
788     return visibleChainAnnotation;
789   }
790
791   public void setVisibleChainAnnotation(boolean visibleChainAnnotation)
792   {
793     this.visibleChainAnnotation = visibleChainAnnotation;
794   }
795
796 }