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