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.xray.parsers;
39
40 import ffx.crystal.Crystal;
41
42 import static ffx.utilities.TinkerUtils.version;
43 import static java.lang.String.format;
44 import static org.apache.commons.math3.util.FastMath.pow;
45 import static org.apache.commons.math3.util.FastMath.sqrt;
46
47 import java.io.DataOutputStream;
48 import java.io.File;
49 import java.io.FileOutputStream;
50 import java.nio.ByteBuffer;
51 import java.nio.ByteOrder;
52 import java.util.logging.Level;
53 import java.util.logging.Logger;
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69 public class CCP4MapWriter {
70
71 private static final Logger logger = Logger.getLogger(CCP4MapWriter.class.getName());
72 private final String filename;
73 private final Crystal crystal;
74 private final int extx, exty, extz;
75 private final int orix, oriy, oriz, nx, ny, nz;
76 private int stride;
77
78
79
80
81
82
83
84
85
86
87 public CCP4MapWriter(int extx, int exty, int extz, Crystal crystal, String filename) {
88 this.orix = 0;
89 this.oriy = 0;
90 this.oriz = 0;
91 this.extx = extx;
92 this.exty = exty;
93 this.extz = extz;
94 this.nx = extx;
95 this.ny = exty;
96 this.nz = extz;
97 this.crystal = crystal;
98 this.filename = filename;
99 this.stride = 2;
100 }
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117 public CCP4MapWriter(
118 int orix, int oriy, int oriz,
119 int extx, int exty, int extz,
120 int nx, int ny, int nz,
121 Crystal crystal, String filename) {
122 this.orix = orix;
123 this.oriy = oriy;
124 this.oriz = oriz;
125 this.extx = extx;
126 this.exty = exty;
127 this.extz = extz;
128 this.nx = nx;
129 this.ny = ny;
130 this.nz = nz;
131 this.crystal = crystal;
132 this.filename = filename;
133 this.stride = 2;
134 }
135
136
137
138
139
140
141 public void setStride(int stride) {
142 this.stride = stride;
143 }
144
145
146
147
148
149
150 public void write(double[] data) {
151 write(data, false);
152 }
153
154
155
156
157
158
159
160 public void write(double[] data, boolean norm) {
161 ByteOrder b = ByteOrder.nativeOrder();
162 FileOutputStream fos;
163 DataOutputStream dos;
164
165 double min = Double.POSITIVE_INFINITY;
166 double max = Double.NEGATIVE_INFINITY;
167 double mean = 0.0;
168 double sd = 0.0;
169
170 int n = 0;
171 for (int k = 0; k < extz; k++) {
172 for (int j = 0; j < exty; j++) {
173 for (int i = 0; i < extx; i++) {
174 int index = stride * (i + extx * (j + exty * k));
175 n++;
176 if (data[index] < min) {
177 min = data[index];
178 }
179 if (data[index] > max) {
180 max = data[index];
181 }
182 mean += (data[index] - mean) / n;
183 }
184 }
185 }
186
187 n = 0;
188 for (int k = 0; k < extz; k++) {
189 for (int j = 0; j < exty; j++) {
190 for (int i = 0; i < extx; i++) {
191 int index = stride * (i + extx * (j + exty * k));
192 sd += pow(data[index] - mean, 2.0);
193 n++;
194 }
195 }
196 }
197 sd = sqrt(sd / n);
198
199 if (norm) {
200 for (int k = 0; k < extz; k++) {
201 for (int j = 0; j < exty; j++) {
202 for (int i = 0; i < extx; i++) {
203 int index = stride * (i + extx * (j + exty * k));
204 data[index] = (data[index] - mean) / sd;
205 }
206 }
207 }
208
209 write(data, false);
210 }
211
212 File file = version(new File(filename));
213 String name = file.getName();
214
215 try {
216 if (logger.isLoggable(Level.INFO)) {
217 StringBuilder sb = new StringBuilder();
218 sb.append(format("\n Writing CCP4 map file: \"%s\"", name));
219 sb.append(format("\n Map min: %g max: %g mean: %g standard dev.: %g", min, max, mean, sd));
220 logger.info(sb.toString());
221 }
222
223 fos = new FileOutputStream(file);
224 dos = new DataOutputStream(fos);
225
226 byte[] bytes = new byte[2048];
227 int offset = 0;
228
229 int imapdata;
230 float fmapdata;
231 String mapstr;
232
233
234 ByteBuffer bb = ByteBuffer.wrap(bytes);
235 bb.order(b).putInt(extx);
236 bb.order(b).putInt(exty);
237 bb.order(b).putInt(extz);
238
239
240 bb.order(b).putInt(2);
241
242 bb.order(b).putInt(orix);
243 bb.order(b).putInt(oriy);
244 bb.order(b).putInt(oriz);
245 bb.order(b).putInt(nx);
246 bb.order(b).putInt(ny);
247 bb.order(b).putInt(nz);
248
249 bb.order(b).putFloat((float) crystal.a);
250 bb.order(b).putFloat((float) crystal.b);
251 bb.order(b).putFloat((float) crystal.c);
252 bb.order(b).putFloat((float) crystal.alpha);
253 bb.order(b).putFloat((float) crystal.beta);
254 bb.order(b).putFloat((float) crystal.gamma);
255
256 bb.order(b).putInt(1);
257 bb.order(b).putInt(2);
258 bb.order(b).putInt(3);
259
260 bb.order(b).putFloat((float) min);
261 bb.order(b).putFloat((float) max);
262 bb.order(b).putFloat((float) mean);
263
264 bb.order(b).putInt(crystal.spaceGroup.number);
265
266
267
268
269 bb.order(b).putInt(80);
270
271 bb.order(b).putInt(0);
272
273 for (int i = 0; i < 12; i++) {
274 bb.order(b).putFloat(0.0f);
275 }
276
277 for (int i = 0; i < 15; i++) {
278 bb.order(b).putInt(0);
279 }
280 dos.write(bytes, offset, 208);
281 bb.rewind();
282
283 mapstr = "MAP ";
284 dos.writeBytes(mapstr);
285
286
287
288 if (ByteOrder.nativeOrder().equals(ByteOrder.LITTLE_ENDIAN)) {
289 imapdata = 0x4441;
290 } else {
291 imapdata = 0x1111;
292 }
293 bb.order(b).putInt(imapdata);
294
295 bb.order(b).putFloat((float) sd);
296
297 bb.order(b).putInt(1);
298 dos.write(bytes, offset, 12);
299
300 StringBuilder sb = new StringBuilder();
301 sb.append("map data from ffx");
302 while (sb.length() < 80) {
303 sb.append(" ");
304 }
305 dos.writeBytes(sb.toString());
306
307 sb = new StringBuilder();
308 while (sb.length() < 80) {
309 sb.append(" ");
310 }
311 for (int i = 0; i < 9; i++) {
312 dos.writeBytes(sb.toString());
313 }
314
315 sb = new StringBuilder();
316 sb.append("x,y,z");
317 while (sb.length() < 80) {
318 sb.append(" ");
319 }
320 dos.writeBytes(sb.toString());
321
322 bb.rewind();
323 for (int k = 0; k < extz; k++) {
324 for (int j = 0; j < exty; j++) {
325 for (int i = 0; i < extx; i++) {
326 int index = stride * (i + extx * (j + exty * k));
327 fmapdata = (float) data[index];
328 bb.order(b).putFloat(fmapdata);
329 if (!bb.hasRemaining()) {
330 dos.write(bytes);
331 bb.rewind();
332 }
333 }
334 }
335 }
336 if (bb.position() > 0) {
337 dos.write(bytes);
338 bb.rewind();
339 }
340
341 dos.close();
342 } catch (Exception e) {
343 String message = "Fatal exception evaluating structure factors.\n";
344 logger.log(Level.SEVERE, message, e);
345 }
346 }
347 }