1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38 package ffx.algorithms.commands;
39
40 import ffx.algorithms.cli.AlgorithmsCommand;
41 import ffx.algorithms.optimize.TitrationManyBody;
42 import ffx.numerics.Potential;
43 import ffx.potential.MolecularAssembly;
44 import ffx.potential.bonded.Atom;
45 import ffx.potential.bonded.Polymer;
46 import ffx.potential.bonded.Residue;
47 import ffx.potential.bonded.Rotamer;
48 import ffx.potential.bonded.RotamerLibrary;
49 import ffx.potential.bonded.RotamerLibrary.NucleicSugarPucker;
50 import ffx.potential.parameters.TitrationUtils;
51 import ffx.potential.parsers.PDBFilter;
52 import ffx.utilities.FFXBinding;
53 import org.apache.commons.configuration2.CompositeConfiguration;
54 import picocli.CommandLine.Command;
55 import picocli.CommandLine.Option;
56 import picocli.CommandLine.Parameters;
57
58 import java.io.File;
59 import java.util.ArrayList;
60 import java.util.Collections;
61 import java.util.HashSet;
62 import java.util.List;
63 import java.util.Set;
64
65 import static java.lang.String.format;
66 import static org.apache.commons.io.FilenameUtils.getExtension;
67 import static org.apache.commons.io.FilenameUtils.removeExtension;
68
69
70
71
72
73
74
75
76 @Command(description = " Save out rotamers.", name = "SaveRotamers")
77 public class SaveRotamers extends AlgorithmsCommand {
78
79 @Option(names = {"--chain", "-c"}, paramLabel = " ",
80 description = "Single character chain name.")
81 private char c = ' ';
82
83 @Option(names = {"--library", "-l"}, paramLabel = "1", defaultValue = "1",
84 description = "Available rotamer libraries are (1) Ponder and Richards or (2) Richardson.")
85 private int library;
86
87 @Option(names = {"--resid", "-r"}, paramLabel = "1", defaultValue = "1",
88 description = "Residue number.")
89 private int resID;
90
91 @Option(names = {"--independent", "-i"}, paramLabel = "false", defaultValue = "false",
92 description = "Independent draws nucleic acid rotamers independently of chain context.")
93 private boolean independent;
94
95 @Option(names = {"--start", "-s"}, paramLabel = "0", defaultValue = "0",
96 description = "First rotamer to draw (indexed from rotamer 0).")
97 private int start;
98
99 @Option(names = {"--finish", "-f"}, paramLabel = "-1", defaultValue = "-1",
100 description = "Last rotamer to draw (indexed from rotamer 0).")
101 private int finish;
102
103 @Option(names = {"--all", "-x"}, paramLabel = "-1", defaultValue = "-1",
104 description = "Draw all rotamers beginning from the passed index (overrides other options).")
105 private int all;
106
107 @Option(names = {"--upstreamPucker", "-u"}, paramLabel = "false", defaultValue = "false",
108 description = "Adjusts the pucker of the 5' residue to match the rotamer.")
109 private boolean upstreamPucker;
110
111 @Option(names = {"--tR", "--titrateResidue"}, paramLabel = "false", defaultValue = "false",
112 description = "Titrate residues.")
113 private boolean titrateResidue;
114
115 @Parameters(arity = "1", paramLabel = "file",
116 description = "The atomic coordinate file in XYZ or PDB format.")
117 private String filename = null;
118
119 private TitrationManyBody titrationManyBody;
120 private TitrationUtils titrationUtils;
121
122 public SaveRotamers() {
123 super();
124 }
125
126 public SaveRotamers(FFXBinding binding) {
127 super(binding);
128 }
129
130 public SaveRotamers(String[] args) {
131 super(args);
132 }
133
134 @Override
135 public SaveRotamers run() {
136 if (!init()) {
137 return this;
138 }
139
140 activeAssembly = getActiveAssembly(filename);
141 if (activeAssembly == null) {
142 logger.info(helpString());
143 return this;
144 }
145
146 filename = activeAssembly.getFile().getAbsolutePath();
147 CompositeConfiguration properties = activeAssembly.getProperties();
148
149 if (titrateResidue) {
150 List<Residue> residues = activeAssembly.getResidueList();
151 List<Integer> resNumberList = new ArrayList<>();
152 for (Residue residue : residues) {
153 resNumberList.add(residue.getResidueNumber());
154 }
155 titrationManyBody = new TitrationManyBody(filename, activeAssembly.getForceField(),
156 resNumberList, 7.0);
157 MolecularAssembly protonatedAssembly = titrationManyBody.getProtonatedAssembly();
158 setActiveAssembly(protonatedAssembly);
159 }
160 RotamerLibrary rLib = new RotamerLibrary(
161 RotamerLibrary.ProteinLibrary.intToProteinLibrary(library), true);
162 String chain = String.valueOf(c);
163
164 boolean saveAllRotamers = false;
165 int allStart = 0;
166 if (all > -1) {
167 allStart = all;
168 saveAllRotamers = true;
169 }
170
171 logger.info("\n Saving rotamers for residue number " + resID + " of chain " + chain + ".");
172
173 RotamerLibrary.initializeDefaultAtomicCoordinates(activeAssembly.getChains());
174 Polymer polymer = activeAssembly.getChain(chain);
175 if (polymer == null) {
176 logger.info(" Polymer + " + chain + " does not exist.");
177 return this;
178 }
179 Residue residue = polymer.getResidue(resID);
180 if (residue == null) {
181 logger.info(" Residue + " + resID + " does not exist.");
182 return this;
183 }
184
185 residue.setRotamers(rLib);
186 Rotamer[] rotamers = residue.getRotamers();
187 if (rotamers == null) {
188 logger.severe(" There are no rotamers for residue + " + residue.toString());
189 }
190 int nrotamers = rotamers.length;
191
192 boolean isDeoxy = false;
193 Residue prevResidue = null;
194 if (upstreamPucker) {
195 try {
196 if (residue.getResidueType() == Residue.ResidueType.NA) {
197 prevResidue = residue.getPreviousResidue();
198 if (prevResidue == null) {
199 upstreamPucker = false;
200 } else {
201 Atom HOs = (Atom) prevResidue.getAtomNode("HO'");
202 if (HOs == null) {
203 isDeoxy = true;
204 }
205 }
206 } else {
207 upstreamPucker = false;
208 }
209 } catch (Exception e) {
210 upstreamPucker = false;
211 }
212 }
213
214 String ext = getExtension(filename);
215 filename = removeExtension(filename);
216 Set<Atom> excludeAtoms = new HashSet<>();
217 List<Residue> removeAtomList = new ArrayList<>();
218 if (saveAllRotamers) {
219 if (allStart >= nrotamers) {
220 logger.info(" Specified start range is outside of rotamer range. No action taken.");
221 } else {
222 for (int i = allStart; i < nrotamers; i++) {
223 RotamerLibrary.applyRotamer(residue, rotamers[i], independent);
224 logger.info(rotamers[i].getName());
225 if (upstreamPucker) {
226 double prevDelta = rotamers[i].chi1;
227 NucleicSugarPucker sugarPucker = NucleicSugarPucker.checkPucker(prevDelta, isDeoxy);
228 RotamerLibrary.applySugarPucker(prevResidue, sugarPucker, isDeoxy, true);
229 }
230 if (ext.toUpperCase().contains("XYZ")) {
231 logger.info(format(" Saving rotamer %d", i));
232 algorithmFunctions.saveAsXYZ(activeAssembly, new File(filename + ".xyz"));
233 } else {
234 if (titrateResidue) {
235 excludeAtoms.clear();
236 removeAtomList.clear();
237 removeAtomList.add(residue);
238
239 int[] currentRot = new int[]{i};
240 titrationManyBody.excludeExcessAtoms(excludeAtoms, currentRot, removeAtomList);
241 File modelFile = new File(filename + ".pdb");
242 PDBFilter pdbFilter = new PDBFilter(modelFile, activeAssembly, activeAssembly.getForceField(),
243 properties);
244 if (!pdbFilter.writeFile(modelFile, false, excludeAtoms, true, true)) {
245 logger.info(format(" Save failed for %s", activeAssembly));
246 }
247 } else {
248 logger.info(format(" Saving rotamer %d", i));
249 algorithmFunctions.saveAsPDB(activeAssembly, new File(filename + ".pdb"));
250 }
251 }
252 }
253 }
254 } else {
255 if (start >= nrotamers) {
256 logger.info(" Specified start range is outside of rotamer range. No action taken.");
257 } else {
258 if (finish >= nrotamers) {
259 finish = nrotamers - 1;
260 } else if (finish < start) {
261 logger.info(
262 " Specified finish point is before the start point drawing only rotamer " + start);
263 finish = start;
264 }
265 for (int i = start; i <= finish; i++) {
266 RotamerLibrary.applyRotamer(residue, rotamers[i], independent);
267 if (upstreamPucker) {
268 double prevDelta = rotamers[i].chi1;
269 NucleicSugarPucker sugarPucker = NucleicSugarPucker.checkPucker(prevDelta, isDeoxy);
270 RotamerLibrary.applySugarPucker(prevResidue, sugarPucker, isDeoxy, true);
271 }
272 if (ext.toUpperCase().contains("XYZ")) {
273 logger.info(format(" Saving rotamer %d", i));
274 algorithmFunctions.saveAsXYZ(activeAssembly, new File(filename + ".xyz"));
275 } else {
276 logger.info(format(" Saving rotamer %d", i));
277 algorithmFunctions.saveAsPDB(activeAssembly, new File(filename + ".pdb"));
278 }
279 }
280 }
281 }
282
283 return this;
284 }
285
286 @Override
287 public List<Potential> getPotentials() {
288 return Collections.emptyList();
289 }
290 }