render_hkl.c 19.6 KB
Newer Older
1
2
3
4
5
/*
 * render_hkl.c
 *
 * Draw pretty renderings of reflection lists
 *
6
7
8
9
 * Copyright © 2012 Deutsches Elektronen-Synchrotron DESY,
 *                  a research centre of the Helmholtz Association.
 *
 * Authors:
10
 *   2010-2012 Thomas White <taw@physics.org>
11
 *
Thomas White's avatar
Thomas White committed
12
13
14
15
16
17
18
19
20
21
22
23
24
25
 * This file is part of CrystFEL.
 *
 * CrystFEL is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * CrystFEL is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with CrystFEL.  If not, see <http://www.gnu.org/licenses/>.
26
27
28
29
30
31
32
33
34
35
36
37
38
39
 *
 */


#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

#include <stdarg.h>
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <unistd.h>
#include <getopt.h>
40
#ifdef HAVE_CAIRO
Thomas White's avatar
Thomas White committed
41
42
#include <cairo.h>
#include <cairo-pdf.h>
43
#endif
44
45

#include "utils.h"
46
#include "symmetry.h"
47
#include "render.h"
48
#include "render_hkl.h"
Thomas White's avatar
Thomas White committed
49
50
#include "reflist.h"
#include "reflist-utils.h"
51
52


53
54
55
#define KEY_FILENAME "key.pdf"


56
57
58
59
static void show_help(const char *s)
{
	printf("Syntax: %s [options] <file.hkl>\n\n", s);
	printf(
60
"Render intensity lists in 2D slices.\n"
61
"\n"
62
"  -d, --down=<h>,<k>,<l>  Indices for the axis in the downward direction.\n"
63
"                           Default: 1,0,0.\n"
64
"  -r, --right=<h>,<k>,<l> Indices for the axis in the 'right' (roughly)\n"
65
"                           direction.  Default: 0,1,0.\n"
Thomas White's avatar
Thomas White committed
66
"  -o, --output=<filename> Output filename.  Default: za.pdf\n"
67
68
69
"      --boost=<val>       Squash colour scale by <val>.\n"
"  -p, --pdb=<file>        PDB file from which to get the unit cell.\n"
"  -y, --symmetry=<sym>    Expand reflections according to point group <sym>.\n"
70
71
72
73
"\n"
"  -c, --colscale=<scale>  Use the given colour scale.  Choose from:\n"
"                           mono    : Greyscale, black is zero.\n"
"                           invmono : Greyscale, white is zero.\n"
Thomas White's avatar
Thomas White committed
74
"                           colour  : Colour scale:\n"
75
76
"                                     black-blue-pink-red-orange-yellow-white\n"
"\n"
77
78
"  -w  --weighting=<wght>  Colour/shade the reciprocal lattice points\n"
"                           according to:\n"
79
80
"                            I      : the intensity of the reflection.\n"
"                            sqrtI  : the square root of the intensity.\n"
Thomas White's avatar
Thomas White committed
81
"                            count  : the number of measurements for the reflection.\n"
82
"                                     (after correcting for 'epsilon')\n"
Thomas White's avatar
Thomas White committed
83
"                            rawcts : the raw number of measurements for the\n"
84
"                                     reflection (no 'epsilon' correction).\n"
85
86
"\n"
"      --colour-key        Draw (only) the key for the current colour scale.\n"
Thomas White's avatar
Thomas White committed
87
88
"                           The key will be written to 'key.pdf' in the\n"
"                           current directory.\n"
Thomas White's avatar
Thomas White committed
89
"\n"
90
"  -h, --help              Display this help message.\n"
91
);
92
93
94
}


95
#ifdef HAVE_CAIRO
96
97


98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
static double max_value(RefList *list, int wght, const SymOpList *sym)
{
	Reflection *refl;
	RefListIterator *iter;
	double max = -INFINITY;
	SymOpMask *m;

	m = new_symopmask(sym);

	for ( refl = first_refl(list, &iter);
	      refl != NULL;
	      refl = next_refl(refl, iter) )
	{
		double val;
		int n;
		signed int h, k, l;

		get_indices(refl, &h, &k, &l);

		special_position(sym, m, h, k, l);
		n = num_equivs(sym, m);

Thomas White's avatar
Thomas White committed
120
		switch ( wght ) {
Thomas White's avatar
Thomas White committed
121
122

			case WGHT_I :
123
124
			val = get_intensity(refl);
			break;
Thomas White's avatar
Thomas White committed
125
126

			case WGHT_SQRTI :
127
128
129
			val = get_intensity(refl);
			val = (val>0.0) ? sqrt(val) : 0.0;
			break;
Thomas White's avatar
Thomas White committed
130
131

			case WGHT_COUNTS :
132
133
134
			val = get_redundancy(refl);
			val /= (double)n;
			break;
Thomas White's avatar
Thomas White committed
135
136

			case WGHT_RAWCOUNTS :
137
138
			val = get_redundancy(refl);
			break;
Thomas White's avatar
Thomas White committed
139
140

			default :
141
142
			ERROR("Invalid weighting.\n");
			abort();
Thomas White's avatar
Thomas White committed
143

144
145
146
147
148
149
150
151
152
		}

		if ( val > max ) max = val;
	}

	return max;
}


Thomas White's avatar
Thomas White committed
153
154
static void draw_circles(double xh, double xk, double xl,
                         double yh, double yk, double yl,
155
                         signed int zh, signed int zk, signed int zl,
Thomas White's avatar
Thomas White committed
156
                         RefList *list, const SymOpList *sym,
157
158
159
                         cairo_t *dctx, int wght, double boost, int colscale,
                         UnitCell *cell, double radius, double theta,
                         double as, double bs, double cx, double cy,
160
                         double scale, double max_val)
161
{
Thomas White's avatar
Thomas White committed
162
163
	Reflection *refl;
	RefListIterator *iter;
Thomas White's avatar
Thomas White committed
164
	SymOpMask *m;
165
	double nx, ny;
166

Thomas White's avatar
Thomas White committed
167
168
	m = new_symopmask(sym);

169
170
171
	nx = xh*xh + xk*xk + xl*xl;
	ny = yh*yh + yk*yk + yl*yl;

Thomas White's avatar
Thomas White committed
172
173
174
	/* Iterate over all reflections */
	for ( refl = first_refl(list, &iter);
	      refl != NULL;
175
176
	      refl = next_refl(refl, iter) )
	{
Thomas White's avatar
Thomas White committed
177
		double u, v, val;
Thomas White's avatar
Thomas White committed
178
		signed int ha, ka, la;
Thomas White's avatar
Thomas White committed
179
		int i, n;
180
		double r, g, b;
Thomas White's avatar
Thomas White committed
181
182
183

		get_indices(refl, &ha, &ka, &la);

Thomas White's avatar
Thomas White committed
184
185
		special_position(sym, m, ha, ka, la);
		n = num_equivs(sym, m);
Thomas White's avatar
Thomas White committed
186
187

		for ( i=0; i<n; i++ ) {
Thomas White's avatar
Thomas White committed
188
189

			signed int h, k, l;
Thomas White's avatar
Thomas White committed
190
			double xi, yi;
Thomas White's avatar
Thomas White committed
191

Thomas White's avatar
Thomas White committed
192
			get_equiv(sym, m, i, ha, ka, la, &h, &k, &l);
Thomas White's avatar
Thomas White committed
193
194
195
196

			/* Is the reflection in the zone? */
			if ( h*zh + k*zk + l*zl != 0 ) continue;

197
198
			xi = (h*xh + k*xk + l*xl) / nx;
			yi = (h*yh + k*yk + l*yl) / ny;
Thomas White's avatar
Thomas White committed
199
200

			switch ( wght) {
Thomas White's avatar
Thomas White committed
201
202

				case WGHT_I :
Thomas White's avatar
Thomas White committed
203
204
				val = get_intensity(refl);
				break;
Thomas White's avatar
Thomas White committed
205
206

				case WGHT_SQRTI :
Thomas White's avatar
Thomas White committed
207
208
209
				val = get_intensity(refl);
				val = (val>0.0) ? sqrt(val) : 0.0;
				break;
Thomas White's avatar
Thomas White committed
210
211

				case WGHT_COUNTS :
Thomas White's avatar
Thomas White committed
212
				val = get_redundancy(refl);
Thomas White's avatar
Thomas White committed
213
				val /= (double)n;
Thomas White's avatar
Thomas White committed
214
				break;
Thomas White's avatar
Thomas White committed
215
216

				case WGHT_RAWCOUNTS :
Thomas White's avatar
Thomas White committed
217
218
				val = get_redundancy(refl);
				break;
Thomas White's avatar
Thomas White committed
219
220

				default :
Thomas White's avatar
Thomas White committed
221
222
				ERROR("Invalid weighting.\n");
				abort();
Thomas White's avatar
Thomas White committed
223

Thomas White's avatar
Thomas White committed
224
			}
225

Thomas White's avatar
Thomas White committed
226
227
228
			/* Absolute location in image based on 2D basis */
			u = (double)xi*as*sin(theta);
			v = (double)xi*as*cos(theta) + (double)yi*bs;
229

230
231
232
			cairo_arc(dctx, ((double)cx)+u*scale,
				        ((double)cy)+v*scale,
				        radius, 0.0, 2.0*M_PI);
233

234
235
236
237
			render_scale(val, max_val/boost, colscale,
				     &r, &g, &b);
			cairo_set_source_rgb(dctx, r, g, b);
			cairo_fill(dctx);
Thomas White's avatar
Thomas White committed
238

239
240
241
		}

	}
Thomas White's avatar
Thomas White committed
242
243

	free_symopmask(m);
244
245
246
}


247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
static void render_overlined_indices(cairo_t *dctx,
                                     signed int h, signed int k, signed int l)
{
	char tmp[256];
	cairo_text_extents_t size;
	double x, y;
	const double sh = 39.0;

	cairo_get_current_point(dctx, &x, &y);
	cairo_set_line_width(dctx, 4.0);

	/* Draw 'h' */
	snprintf(tmp, 255, "%i", abs(h));
	cairo_text_extents(dctx, tmp, &size);
	cairo_show_text(dctx, tmp);
	cairo_fill(dctx);
	if ( h < 0 ) {
		cairo_move_to(dctx, x+size.x_bearing, y-sh);
		cairo_rel_line_to(dctx, size.width, 0.0);
		cairo_stroke(dctx);
	}
	x += size.x_advance;

	/* Draw 'k' */
	cairo_move_to(dctx, x, y);
	snprintf(tmp, 255, "%i", abs(k));
	cairo_text_extents(dctx, tmp, &size);
	cairo_show_text(dctx, tmp);
	cairo_fill(dctx);
	if ( k < 0 ) {
		cairo_move_to(dctx, x+size.x_bearing, y-sh);
		cairo_rel_line_to(dctx, size.width, 0.0);
		cairo_stroke(dctx);
	}
	x += size.x_advance;

	/* Draw 'l' */
	cairo_move_to(dctx, x, y);
	snprintf(tmp, 255, "%i", abs(l));
	cairo_text_extents(dctx, tmp, &size);
	cairo_show_text(dctx, tmp);
	cairo_fill(dctx);
	if ( l < 0 ) {
		cairo_move_to(dctx, x+size.x_bearing, y-sh);
		cairo_rel_line_to(dctx, size.width, 0.0);
		cairo_stroke(dctx);
	}
}


Thomas White's avatar
Thomas White committed
297
static void render_za(UnitCell *cell, RefList *list,
Thomas White's avatar
Thomas White committed
298
299
                      double boost, const SymOpList *sym, int wght,
                      int colscale,
300
                      signed int xh, signed int xk, signed int xl,
301
                      signed int yh, signed int yk, signed int yl,
302
                      const char *outfile, double scale_top)
303
{
Thomas White's avatar
Thomas White committed
304
305
	cairo_surface_t *surface;
	cairo_t *dctx;
Thomas White's avatar
Thomas White committed
306
307
	double max_val;
	double scale1, scale2, scale;
Thomas White's avatar
Thomas White committed
308
	double sep_u, sep_v, max_r;
309
	double u, v;
Thomas White's avatar
Thomas White committed
310
	double as, bs, theta;
311
312
313
	double asx, asy, asz;
	double bsx, bsy, bsz;
	double csx, csy, csz;
Thomas White's avatar
Thomas White committed
314
	float wh, ht;
315
316
317
318
319
320
321
	signed int zh, zk, zl;
	double xx, xy, xz;
	double yx, yy, yz;
	char tmp[256];
	cairo_text_extents_t size;
	double cx, cy;
	const double border = 200.0;
322
	int png;
Thomas White's avatar
Thomas White committed
323
	double rmin, rmax;
324

325
	/* Vector product to determine the zone axis. */
326
327
328
329
	zh = xk*yl - xl*yk;
	zk = - xh*yl + xl*yh;
	zl = xh*yk - xk*yh;
	STATUS("Zone axis is %i %i %i\n", zh, zk, zl);
330

331
	/* Size of output and centre definition */
Thomas White's avatar
Thomas White committed
332
333
	wh = 1024;
	ht = 1024;
334

Thomas White's avatar
Thomas White committed
335
	/* Work out reciprocal lattice spacings and angles for this cut */
336
	if ( cell_get_reciprocal(cell, &asx, &asy, &asz,
Thomas White's avatar
Thomas White committed
337
	                          &bsx, &bsy, &bsz,
338
339
340
341
	                          &csx, &csy, &csz) ) {
		ERROR("Couldn't get reciprocal parameters\n");
		return;
	}
342
343
344
345
346
347
348
	xx = xh*asx + xk*bsx + xl*csx;
	xy = xh*asy + xk*bsy + xl*csy;
	xz = xh*asz + xk*bsz + xl*csz;
	yx = yh*asx + yk*bsx + yl*csx;
	yy = yh*asy + yk*bsy + yl*csy;
	yz = yh*asz + yk*bsz + yl*csz;
	theta = angle_between(xx, xy, xz, yx, yy, yz);
Thomas White's avatar
Thomas White committed
349
350
	as = modulus(xx, xy, xz);
	bs = modulus(yx, yy, yz);
351

Thomas White's avatar
Thomas White committed
352
353
354
355
	resolution_limits(list, cell, &rmin, &rmax);
	printf("Resolution limits: 1/d = %.2f - %.2f nm^-1"
	       " (d = %.2f - %.2f A)\n",
	       rmin/1e9, rmax/1e9, (1.0/rmin)/1e-10, (1.0/rmax)/1e-10);
356

357
	max_val = max_value(list, wght, sym);
358
	if ( max_val <= 0.0 ) {
359
		STATUS("Couldn't find max value.\n");
Thomas White's avatar
Thomas White committed
360
		return;
361
362
	}

363
364
365
366
367
	/* Use manual scale top if specified */
	if ( scale_top > 0.0 ) {
		max_val = scale_top;
	}

Thomas White's avatar
Thomas White committed
368
369
370
	scale1 = ((double)wh-border) / (2.0*rmax);
	scale2 = ((double)ht-border) / (2.0*rmax);
	scale = (scale1 < scale2) ? scale1 : scale2;
371

Thomas White's avatar
Thomas White committed
372
	/* Work out the spot radius */
373
374
	sep_u = scale*as;
	sep_v = scale*bs;
375
	max_r = (sep_u < sep_v) ? sep_u : sep_v;
376
	max_r /= 2.0;  /* Max radius is half the separation */
377
	max_r -= (max_r/10.0);  /* Add a tiny separation between circles */
378

Thomas White's avatar
Thomas White committed
379
	/* Create surface */
380
381
382
383
384
385
386
387
	if ( strcmp(outfile+strlen(outfile)-4, ".png") == 0 ) {
		png = 1;
		surface = cairo_image_surface_create(CAIRO_FORMAT_ARGB32,
		                                     wh, ht);
	} else {
		png = 0;
		surface = cairo_pdf_surface_create(outfile, wh, ht);
	}
388
389

	if ( cairo_surface_status(surface) != CAIRO_STATUS_SUCCESS ) {
Thomas White's avatar
Thomas White committed
390
		ERROR("Couldn't create Cairo surface\n");
391
392
393
394
395
		cairo_surface_destroy(surface);
		return;
	}

	dctx = cairo_create(surface);
Thomas White's avatar
Thomas White committed
396
397
398
399
400
	if ( cairo_status(dctx) != CAIRO_STATUS_SUCCESS ) {
		ERROR("Couldn't create Cairo context\n");
		cairo_surface_destroy(surface);
		return;
	}
401
402
403
404
405
406
407
408
409
410
411
412
413
414

	/* Black background */
	cairo_rectangle(dctx, 0.0, 0.0, wh, ht);
	cairo_set_source_rgb(dctx, 0.0, 0.0, 0.0);
	cairo_fill(dctx);

	/* Test size of text that goes to the right(ish) */
	cairo_set_font_size(dctx, 40.0);
	snprintf(tmp, 255, "%i%i%i", abs(xh), abs(xk), abs(xl));
	cairo_text_extents(dctx, tmp, &size);

	cx = 532.0 - size.width;
	cy = 512.0 - 20.0;

415
	draw_circles(xh, xk, xl, yh, yk, yl, zh, zk, zl,
Thomas White's avatar
Thomas White committed
416
	             list, sym, dctx, wght, boost, colscale, cell,
417
	             max_r, theta, as, bs, cx, cy, scale,
418
	             max_val);
419

Thomas White's avatar
Thomas White committed
420
	/* Centre marker */
421
422
	cairo_arc(dctx, (double)cx,
			(double)cy, max_r, 0, 2*M_PI);
Thomas White's avatar
Thomas White committed
423
424
425
	cairo_set_source_rgb(dctx, 1.0, 0.0, 0.0);
	cairo_fill(dctx);

426
	/* Draw indexing lines */
427
	cairo_set_line_cap(dctx, CAIRO_LINE_CAP_ROUND);
428
	cairo_set_line_width(dctx, 4.0);
429
	cairo_move_to(dctx, (double)cx, (double)cy);
Thomas White's avatar
Thomas White committed
430
431
	u = rmax*sin(theta);
	v = rmax*cos(theta);
432
	cairo_line_to(dctx, cx+u*scale, cy+v*scale);
433
434
435
	cairo_set_source_rgb(dctx, 0.0, 1.0, 0.0);
	cairo_stroke(dctx);

436
	cairo_set_font_size(dctx, 40.0);
437
	snprintf(tmp, 255, "%i%i%i", abs(xh), abs(xk), abs(xl));
438
439
440
	cairo_text_extents(dctx, tmp, &size);

	cairo_move_to(dctx, cx+u*scale + 20.0, cy+v*scale + size.height/2.0);
441
	render_overlined_indices(dctx, xh, xk, xl);
442
443
	cairo_fill(dctx);

444
	snprintf(tmp, 255, "%i%i%i", abs(yh), abs(yk), abs(yl));
445
446
447
	cairo_text_extents(dctx, tmp, &size);

	cairo_move_to(dctx, (double)cx, (double)cy);
Thomas White's avatar
Thomas White committed
448
	cairo_line_to(dctx, cx, cy+rmax*scale);
449
450
451
	cairo_set_source_rgb(dctx, 0.0, 1.0, 0.0);
	cairo_stroke(dctx);

Thomas White's avatar
Thomas White committed
452
453
	cairo_move_to(dctx, cx - size.width/2.0,
	                    cy+rmax*scale + size.height + 20.0);
454
	render_overlined_indices(dctx, yh, yk, yl);
455
456
	cairo_fill(dctx);

457
458
459
460
461
462
463
	if ( png ) {
		int r = cairo_surface_write_to_png(surface, outfile);
		if ( r != CAIRO_STATUS_SUCCESS ) {
			ERROR("Failed to write PNG to '%s'\n", outfile);
		}
	}

Thomas White's avatar
Thomas White committed
464
465
466
	cairo_surface_finish(surface);
	cairo_destroy(dctx);
}
467

468

469
static int render_key(int colscale, double scale_top)
470
471
472
{
	cairo_surface_t *surface;
	cairo_t *dctx;
473
474
	double top, wh, ht, y;
	double slice;
475

476
477
478
479
480
481
482
483
484
	wh = 128.0;
	ht = 1024.0;
	slice = 1.0;

	if ( scale_top > 0.0 ) {
		top = scale_top;
	} else {
		top = 1.0;
	}
485

486
	surface = cairo_pdf_surface_create(KEY_FILENAME, wh, ht);
487
488
489
490
491
492
493
494
495

	if ( cairo_surface_status(surface) != CAIRO_STATUS_SUCCESS ) {
		fprintf(stderr, "Couldn't create Cairo surface\n");
		cairo_surface_destroy(surface);
		return 1;
	}

	dctx = cairo_create(surface);

496
	for ( y=0.0; y<ht; y+=slice ) {
497

498
499
500
		double r, g, b;
		double val;
		double v = y;
501

502
503
504
505
506
507
508
509
510
511
512
		cairo_rectangle(dctx, 0.0, ht-y, wh/2.0, slice);

		if ( colscale == SCALE_RATIO ) {
			if ( v < ht/2.0 ) {
				val = v/(ht/2.0);
			} else {
				val = (((v-ht/2.0)/(ht/2.0))*(top-1.0))+1.0;
			}
		} else {
			val = v/ht;
		}
513

514
		render_scale(val, top, colscale, &r, &g, &b);
515
516
		cairo_set_source_rgb(dctx, r, g, b);

517
		cairo_stroke_preserve(dctx);
518
519
520
521
		cairo_fill(dctx);

	}

522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
	if ( colscale == SCALE_RATIO ) {

		cairo_text_extents_t size;
		char tmp[32];

		cairo_rectangle(dctx, 0.0, ht/2.0-2.0, wh/2.0, 4.0);
		cairo_set_source_rgb(dctx, 0.0, 0.0, 0.0);
		cairo_stroke_preserve(dctx);
		cairo_fill(dctx);

		cairo_set_font_size(dctx, 20.0);
		cairo_text_extents(dctx, "1.0", &size);
		cairo_move_to(dctx, wh/2.0+5.0, ht/2.0+size.height/2.0);
		cairo_show_text(dctx, "1.0");

		cairo_set_font_size(dctx, 20.0);
		cairo_text_extents(dctx, "0.0", &size);
		cairo_move_to(dctx, wh/2.0+5.0, ht-5.0);
		cairo_show_text(dctx, "0.0");

		cairo_set_font_size(dctx, 20.0);
		snprintf(tmp, 31, "%.1f", top);
		cairo_text_extents(dctx, tmp, &size);
		cairo_move_to(dctx, wh/2.0+5.0, size.height+5.0);
		cairo_show_text(dctx, tmp);

	}


551
552
553
	cairo_surface_finish(surface);
	cairo_destroy(dctx);

554
555
	STATUS("Colour key written to "KEY_FILENAME"\n");

556
557
	return 0;
}
558
559
560
561
562


#else  /* HAVE_CAIRO */


563
static int render_key(int colscale, double scale_top)
564
565
566
567
568
569
570
571
{
	ERROR("This version of CrystFEL was compiled without Cairo");
	ERROR(" support, which is required to draw the colour");
	ERROR(" scale.  Sorry!\n");
	return 1;
}


572
static void render_za(UnitCell *cell, RefList *list,
573
574
575
                      double boost, const char *sym, int wght, int colscale,
                      signed int xh, signed int xk, signed int xl,
                      signed int yh, signed int yk, signed int yl,
576
                      const char *outfile, double scale_top)
577
578
579
580
581
582
583
584
{
	ERROR("This version of CrystFEL was compiled without Cairo");
	ERROR(" support, which is required to plot a zone axis");
	ERROR(" pattern.  Sorry!\n");
}


#endif /* HAVE_CAIRO */
Thomas White's avatar
Thomas White committed
585
586
587
588
589
590


int main(int argc, char *argv[])
{
	int c;
	UnitCell *cell;
Thomas White's avatar
Thomas White committed
591
	RefList *list;
Thomas White's avatar
Thomas White committed
592
	char *infile;
Thomas White's avatar
Thomas White committed
593
	int config_sqrt = 0;
594
	int config_colkey = 0;
595
	int config_zawhinge = 0;
596
	char *pdb = NULL;
597
	int r = 0;
598
	double boost = 1.0;
Thomas White's avatar
Thomas White committed
599
600
	char *sym_str = NULL;
	SymOpList *sym;
601
602
	char *weighting = NULL;
	int wght;
603
604
	int colscale;
	char *cscale = NULL;
605
606
607
608
	signed int dh=1, dk=0, dl=0;
	signed int rh=0, rk=1, rl=0;
	char *down = NULL;
	char *right = NULL;
609
	char *outfile = NULL;
610
611
	double scale_top = -1.0;
	char *endptr;
Thomas White's avatar
Thomas White committed
612
613
614
615

	/* Long options */
	const struct option longopts[] = {
		{"help",               0, NULL,               'h'},
616
		{"zone-axis",          0, &config_zawhinge,    1},
617
		{"output",             1, NULL,               'o'},
618
		{"pdb",                1, NULL,               'p'},
619
		{"boost",              1, NULL,               'b'},
620
		{"symmetry",           1, NULL,               'y'},
621
		{"weighting",          1, NULL,               'w'},
622
		{"colscale",           1, NULL,               'c'},
623
624
		{"down",               1, NULL,               'd'},
		{"right",              1, NULL,               'r'},
625
		{"counts",             0, &config_sqrt,        1},
626
		{"colour-key",         0, &config_colkey,      1},
627
		{"scale-top",          1, NULL,                2},
Thomas White's avatar
Thomas White committed
628
629
630
631
		{0, 0, NULL, 0}
	};

	/* Short options */
Thomas White's avatar
Thomas White committed
632
	while ((c = getopt_long(argc, argv, "hp:w:c:y:d:r:o:",
633
	                        longopts, NULL)) != -1) {
Thomas White's avatar
Thomas White committed
634
635

		switch (c) {
Thomas White's avatar
Thomas White committed
636
637

			case 'h' :
Thomas White's avatar
Thomas White committed
638
639
640
			show_help(argv[0]);
			return 0;

Thomas White's avatar
Thomas White committed
641
			case 'p' :
642
643
644
			pdb = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
645
			case 'b' :
646
647
648
			boost = atof(optarg);
			break;

Thomas White's avatar
Thomas White committed
649
			case 'y' :
Thomas White's avatar
Thomas White committed
650
			sym_str = strdup(optarg);
651
652
			break;

Thomas White's avatar
Thomas White committed
653
			case 'w' :
654
655
656
			weighting = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
657
			case 'c' :
658
659
660
			cscale = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
661
			case 'd' :
662
663
664
			down = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
665
			case 'r' :
666
667
668
			right = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
669
			case 'o' :
670
671
672
			outfile = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
673
			case 2 :
674
675
676
677
678
679
680
681
682
683
684
685
686
687
			errno = 0;
			scale_top = strtod(optarg, &endptr);
			if ( !( (optarg[0] != '\0') && (endptr[0] == '\0') )
			   || (errno != 0) )
			{
				ERROR("Invalid scale top('%s')\n", optarg);
				return 1;
			}
			if ( scale_top < 0.0 ) {
				ERROR("Scale top must be positive.\n");
				return 1;
			}
			break;

Thomas White's avatar
Thomas White committed
688
			case 0 :
Thomas White's avatar
Thomas White committed
689
690
			break;

Thomas White's avatar
Thomas White committed
691
			default :
692
693
			ERROR("Unhandled option '%c'\n", c);
			break;
Thomas White's avatar
Thomas White committed
694

Thomas White's avatar
Thomas White committed
695
696
697
698
		}

	}

699
	if ( config_zawhinge ) {
700
701
702
703
		ERROR("Friendly warning: The --zone-axis option isn't needed"
		      " any longer (I ignored it for you).\n");
	}

704
705
706
	if ( (pdb == NULL) && !config_colkey ) {
		ERROR("You must specify the PDB containing the unit cell.\n");
		return 1;
707
708
	}

Thomas White's avatar
Thomas White committed
709
710
	if ( sym_str == NULL ) {
		sym_str = strdup("1");
711
	}
Thomas White's avatar
Thomas White committed
712
713
	sym = get_pointgroup(sym_str);
	free(sym_str);
714

715
716
717
718
	if ( weighting == NULL ) {
		weighting = strdup("I");
	}

719
	if ( outfile == NULL ) outfile = strdup("za.pdf");
Thomas White's avatar
Thomas White committed
720

721
722
723
724
725
726
	if ( strcmp(weighting, "I") == 0 ) {
		wght = WGHT_I;
	} else if ( strcmp(weighting, "sqrtI") == 0 ) {
		wght = WGHT_SQRTI;
	} else if ( strcmp(weighting, "count") == 0 ) {
		wght = WGHT_COUNTS;
727
728
	} else if ( strcmp(weighting, "counts") == 0 ) {
		wght = WGHT_COUNTS;
729
730
	} else if ( strcmp(weighting, "rawcts") == 0 ) {
		wght = WGHT_RAWCOUNTS;
731
732
733
734
	} else if ( strcmp(weighting, "rawcount") == 0 ) {
		wght = WGHT_RAWCOUNTS;
	} else if ( strcmp(weighting, "rawcounts") == 0 ) {
		wght = WGHT_RAWCOUNTS;
735
736
737
738
	} else {
		ERROR("Unrecognised weighting '%s'\n", weighting);
		return 1;
	}
739
740
741
742
743
744
745
746
747
748
749
750
751
752
	free(weighting);

	if ( cscale == NULL ) {
		cscale = strdup("mono");
	}

	if ( strcmp(cscale, "mono") == 0 ) {
		colscale = SCALE_MONO;
	} else if ( strcmp(cscale, "invmono") == 0 ) {
		colscale = SCALE_INVMONO;
	} else if ( strcmp(cscale, "colour") == 0 ) {
		colscale = SCALE_COLOUR;
	} else if ( strcmp(cscale, "color") == 0 ) {
		colscale = SCALE_COLOUR;
753
754
	} else if ( strcmp(cscale, "ratio") == 0 ) {
		colscale = SCALE_RATIO;
755
756
757
758
759
	} else {
		ERROR("Unrecognised colour scale '%s'\n", cscale);
		return 1;
	}
	free(cscale);
760

761
	if ( config_colkey ) {
762
		return render_key(colscale, scale_top);
763
764
	}

765
766
767
768
769
770
771
772
773
774
775
	if ( (( down == NULL ) && ( right != NULL ))
	  || (( down != NULL ) && ( right == NULL )) ) {
		ERROR("Either specify both 'down' and 'right',"
		      " or neither.\n");
		return 1;
	}
	if ( down != NULL ) {
		int r;
		r = sscanf(down, "%i,%i,%i", &dh, &dk, &dl);
		if ( r != 3 ) {
			ERROR("Invalid format for 'down'\n");
776
777
			return 1;
		}
778
779
780
781
782
783
784
	}
	if ( right != NULL ) {
		int r;
		r = sscanf(right, "%i,%i,%i", &rh, &rk, &rl);
		if ( r != 3 ) {
			ERROR("Invalid format for 'right'\n");
			return 1;
785
786
787
		}
	}

Thomas White's avatar
Thomas White committed
788
789
	infile = argv[optind];

790
	cell = load_cell_from_pdb(pdb);
791
792
793
794
	if ( cell == NULL ) {
		ERROR("Couldn't load unit cell from %s\n", pdb);
		return 1;
	}
Thomas White's avatar
Thomas White committed
795
796
797
	list = read_reflections(infile);
	if ( list == NULL ) {
		ERROR("Couldn't read file '%s'\n", infile);
Thomas White's avatar
Thomas White committed
798
799
		return 1;
	}
Thomas White's avatar
Thomas White committed
800
	if ( check_list_symmetry(list, sym) ) {
Thomas White's avatar
Thomas White committed
801
		ERROR("The input reflection list does not appear to"
Thomas White's avatar
Thomas White committed
802
		      " have symmetry %s\n", symmetry_name(sym));
Thomas White's avatar
Thomas White committed
803
804
		return 1;
	}
Thomas White's avatar
Thomas White committed
805

806
807
	render_za(cell, list, boost, sym, wght, colscale,
	          rh, rk, rl, dh, dk, dl, outfile, scale_top);
Thomas White's avatar
Thomas White committed
808

809
	free(pdb);
Thomas White's avatar
Thomas White committed
810
	free_symoplist(sym);
Thomas White's avatar
Thomas White committed
811
	reflist_free(list);
812
	if ( outfile != NULL ) free(outfile);
813

814
	return r;
815
}