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 activeAssembly = titrationManyBody.getProtonatedAssembly();
158 }
159 RotamerLibrary rLib = new RotamerLibrary(
160 RotamerLibrary.ProteinLibrary.intToProteinLibrary(library), true);
161 String chain = String.valueOf(c);
162
163 boolean saveAllRotamers = false;
164 int allStart = 0;
165 if (all > -1) {
166 allStart = all;
167 saveAllRotamers = true;
168 }
169
170 logger.info("\n Saving rotamers for residue number " + resID + " of chain " + chain + ".");
171
172 RotamerLibrary.initializeDefaultAtomicCoordinates(activeAssembly.getChains());
173 Polymer polymer = activeAssembly.getChain(chain);
174 if (polymer == null) {
175 logger.info(" Polymer + " + chain + " does not exist.");
176 return this;
177 }
178 Residue residue = polymer.getResidue(resID);
179 if (residue == null) {
180 logger.info(" Residue + " + resID + " does not exist.");
181 return this;
182 }
183
184 residue.setRotamers(rLib);
185 Rotamer[] rotamers = residue.getRotamers();
186 if (rotamers == null) {
187 logger.severe(" There are no rotamers for residue + " + residue.toString());
188 }
189 int nrotamers = rotamers.length;
190
191 boolean isDeoxy = false;
192 Residue prevResidue = null;
193 if (upstreamPucker) {
194 try {
195 if (residue.getResidueType() == Residue.ResidueType.NA) {
196 prevResidue = residue.getPreviousResidue();
197 if (prevResidue == null) {
198 upstreamPucker = false;
199 } else {
200 Atom HOs = (Atom) prevResidue.getAtomNode("HO'");
201 if (HOs == null) {
202 isDeoxy = true;
203 }
204 }
205 } else {
206 upstreamPucker = false;
207 }
208 } catch (Exception e) {
209 upstreamPucker = false;
210 }
211 }
212
213 String ext = getExtension(filename);
214 filename = removeExtension(filename);
215 Set<Atom> excludeAtoms = new HashSet<>();
216 List<Residue> removeAtomList = new ArrayList<>();
217 if (saveAllRotamers) {
218 if (allStart >= nrotamers) {
219 logger.info(" Specified start range is outside of rotamer range. No action taken.");
220 } else {
221 for (int i = allStart; i < nrotamers; i++) {
222 RotamerLibrary.applyRotamer(residue, rotamers[i], independent);
223 logger.info(rotamers[i].getName());
224 if (upstreamPucker) {
225 double prevDelta = rotamers[i].chi1;
226 NucleicSugarPucker sugarPucker = NucleicSugarPucker.checkPucker(prevDelta, isDeoxy);
227 RotamerLibrary.applySugarPucker(prevResidue, sugarPucker, isDeoxy, true);
228 }
229 if (ext.toUpperCase().contains("XYZ")) {
230 logger.info(format(" Saving rotamer %d", i));
231 algorithmFunctions.saveAsXYZ(activeAssembly, new File(filename + ".xyz"));
232 } else {
233 if (titrateResidue) {
234 excludeAtoms.clear();
235 removeAtomList.clear();
236 removeAtomList.add(residue);
237
238 int[] currentRot = new int[]{i};
239 titrationManyBody.excludeExcessAtoms(excludeAtoms, currentRot, removeAtomList);
240 File modelFile = new File(filename + ".pdb");
241 PDBFilter pdbFilter = new PDBFilter(modelFile, activeAssembly, activeAssembly.getForceField(),
242 properties);
243 if (!pdbFilter.writeFile(modelFile, false, excludeAtoms, true, true)) {
244 logger.info(format(" Save failed for %s", activeAssembly));
245 }
246 } else {
247 logger.info(format(" Saving rotamer %d", i));
248 algorithmFunctions.saveAsPDB(activeAssembly, new File(filename + ".pdb"));
249 }
250 }
251 }
252 }
253 } else {
254 if (start >= nrotamers) {
255 logger.info(" Specified start range is outside of rotamer range. No action taken.");
256 } else {
257 if (finish >= nrotamers) {
258 finish = nrotamers - 1;
259 } else if (finish < start) {
260 logger.info(
261 " Specified finish point is before the start point drawing only rotamer " + start);
262 finish = start;
263 }
264 for (int i = start; i <= finish; i++) {
265 RotamerLibrary.applyRotamer(residue, rotamers[i], independent);
266 if (upstreamPucker) {
267 double prevDelta = rotamers[i].chi1;
268 NucleicSugarPucker sugarPucker = NucleicSugarPucker.checkPucker(prevDelta, isDeoxy);
269 RotamerLibrary.applySugarPucker(prevResidue, sugarPucker, isDeoxy, true);
270 }
271 if (ext.toUpperCase().contains("XYZ")) {
272 logger.info(format(" Saving rotamer %d", i));
273 algorithmFunctions.saveAsXYZ(activeAssembly, new File(filename + ".xyz"));
274 } else {
275 logger.info(format(" Saving rotamer %d", i));
276 algorithmFunctions.saveAsPDB(activeAssembly, new File(filename + ".pdb"));
277 }
278 }
279 }
280 }
281
282 return this;
283 }
284
285 @Override
286 public List<Potential> getPotentials() {
287 return Collections.emptyList();
288 }
289 }