render_hkl.c 20 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
#include "cell-utils.h"
52
53


54
55
56
#define KEY_FILENAME "key.pdf"


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


96
#ifdef HAVE_CAIRO
97
98


99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
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
121
		switch ( wght ) {
Thomas White's avatar
Thomas White committed
122
123

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

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

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

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

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

145
146
147
148
149
150
151
152
153
		}

		if ( val > max ) max = val;
	}

	return max;
}


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

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

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

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

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

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

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

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

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

			/* Is the reflection in the zone? */
196
			if ( h*zh + k*zk + l*zl != zone ) continue;
Thomas White's avatar
Thomas White committed
197

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

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

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

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

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

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

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

Thomas White's avatar
Thomas White committed
225
			}
226

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

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

235
236
237
238
			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
239

240
241
242
		}

	}
Thomas White's avatar
Thomas White committed
243
244

	free_symopmask(m);
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
297
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
298
static void render_za(UnitCell *cell, RefList *list,
Thomas White's avatar
Thomas White committed
299
300
                      double boost, const SymOpList *sym, int wght,
                      int colscale,
301
                      signed int xh, signed int xk, signed int xl,
302
                      signed int yh, signed int yk, signed int yl,
303
                      const char *outfile, double scale_top, signed int zone)
304
{
Thomas White's avatar
Thomas White committed
305
306
	cairo_surface_t *surface;
	cairo_t *dctx;
Thomas White's avatar
Thomas White committed
307
308
	double max_val;
	double scale1, scale2, scale;
Thomas White's avatar
Thomas White committed
309
	double sep_u, sep_v, max_r;
310
	double u, v;
Thomas White's avatar
Thomas White committed
311
	double as, bs, theta;
312
313
314
	double asx, asy, asz;
	double bsx, bsy, bsz;
	double csx, csy, csz;
Thomas White's avatar
Thomas White committed
315
	float wh, ht;
316
317
318
319
320
321
322
	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;
323
	int png;
Thomas White's avatar
Thomas White committed
324
	double rmin, rmax;
325

326
	/* Vector product to determine the zone axis. */
327
328
329
330
	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);
331

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

Thomas White's avatar
Thomas White committed
336
	/* Work out reciprocal lattice spacings and angles for this cut */
337
	if ( cell_get_reciprocal(cell, &asx, &asy, &asz,
Thomas White's avatar
Thomas White committed
338
	                          &bsx, &bsy, &bsz,
339
340
341
342
	                          &csx, &csy, &csz) ) {
		ERROR("Couldn't get reciprocal parameters\n");
		return;
	}
343
344
345
346
347
348
349
	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
350
351
	as = modulus(xx, xy, xz);
	bs = modulus(yx, yy, yz);
352

Thomas White's avatar
Thomas White committed
353
354
355
356
	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);
357

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

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

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

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

Thomas White's avatar
Thomas White committed
380
	/* Create surface */
381
382
383
384
385
386
387
388
	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);
	}
389
390

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

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

	/* 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;

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

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

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

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

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

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

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

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

458
459
460
461
462
463
464
	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
465
466
467
	cairo_surface_finish(surface);
	cairo_destroy(dctx);
}
468

469

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

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

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

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

	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);

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

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

503
504
505
506
507
508
509
510
511
512
513
		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;
		}
514

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

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

	}

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
551
	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);

	}


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

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

557
558
	return 0;
}
559
560
561
562
563


#else  /* HAVE_CAIRO */


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


573
static void render_za(UnitCell *cell, RefList *list,
574
575
576
                      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,
577
                      const char *outfile, double scale_top)
578
579
580
581
582
583
584
585
{
	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
586
587
588
589
590
591


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

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

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

		switch (c) {
Thomas White's avatar
Thomas White committed
639
640

			case 'h' :
Thomas White's avatar
Thomas White committed
641
642
643
			show_help(argv[0]);
			return 0;

Thomas White's avatar
Thomas White committed
644
			case 'p' :
645
646
647
			pdb = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
648
			case 'b' :
649
650
651
			boost = atof(optarg);
			break;

Thomas White's avatar
Thomas White committed
652
			case 'y' :
Thomas White's avatar
Thomas White committed
653
			sym_str = strdup(optarg);
654
655
			break;

Thomas White's avatar
Thomas White committed
656
			case 'w' :
657
658
659
			weighting = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
660
			case 'c' :
661
662
663
			cscale = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
664
			case 'd' :
665
666
667
			down = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
668
			case 'r' :
669
670
671
			right = strdup(optarg);
			break;

Thomas White's avatar
Thomas White committed
672
			case 'o' :
673
674
675
			outfile = strdup(optarg);
			break;

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

691
692
693
694
695
696
697
698
699
700
701
			case 3 :
			errno = 0;
			zone = strtol(optarg, &endptr, 10);
			if ( !( (optarg[0] != '\0') && (endptr[0] == '\0') )
			   || (errno != 0) )
			{
				ERROR("Invalid zone number ('%s')\n", optarg);
				return 1;
			}
			break;

Thomas White's avatar
Thomas White committed
702
			case 0 :
Thomas White's avatar
Thomas White committed
703
704
			break;

Thomas White's avatar
Thomas White committed
705
			default :
706
707
			ERROR("Unhandled option '%c'\n", c);
			break;
Thomas White's avatar
Thomas White committed
708

Thomas White's avatar
Thomas White committed
709
710
711
712
		}

	}

713
	if ( config_zawhinge ) {
714
715
716
717
		ERROR("Friendly warning: The --zone-axis option isn't needed"
		      " any longer (I ignored it for you).\n");
	}

718
719
720
	if ( (pdb == NULL) && !config_colkey ) {
		ERROR("You must specify the PDB containing the unit cell.\n");
		return 1;
721
722
	}

Thomas White's avatar
Thomas White committed
723
724
	if ( sym_str == NULL ) {
		sym_str = strdup("1");
725
	}
Thomas White's avatar
Thomas White committed
726
727
	sym = get_pointgroup(sym_str);
	free(sym_str);
728

729
730
731
732
	if ( weighting == NULL ) {
		weighting = strdup("I");
	}

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

735
736
737
738
739
740
	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;
741
742
	} else if ( strcmp(weighting, "counts") == 0 ) {
		wght = WGHT_COUNTS;
743
744
	} else if ( strcmp(weighting, "rawcts") == 0 ) {
		wght = WGHT_RAWCOUNTS;
745
746
747
748
	} else if ( strcmp(weighting, "rawcount") == 0 ) {
		wght = WGHT_RAWCOUNTS;
	} else if ( strcmp(weighting, "rawcounts") == 0 ) {
		wght = WGHT_RAWCOUNTS;
749
750
751
752
	} else {
		ERROR("Unrecognised weighting '%s'\n", weighting);
		return 1;
	}
753
754
755
756
757
758
759
760
761
762
763
764
765
766
	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;
767
768
	} else if ( strcmp(cscale, "ratio") == 0 ) {
		colscale = SCALE_RATIO;
769
770
771
772
773
	} else {
		ERROR("Unrecognised colour scale '%s'\n", cscale);
		return 1;
	}
	free(cscale);
774

775
	if ( config_colkey ) {
776
		return render_key(colscale, scale_top);
777
778
	}

779
780
781
782
783
784
785
786
787
788
789
	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");
790
791
			return 1;
		}
792
793
794
795
796
797
798
	}
	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;
799
800
801
		}
	}

Thomas White's avatar
Thomas White committed
802
803
	infile = argv[optind];

804
	cell = load_cell_from_pdb(pdb);
805
806
807
808
	if ( cell == NULL ) {
		ERROR("Couldn't load unit cell from %s\n", pdb);
		return 1;
	}
Thomas White's avatar
Thomas White committed
809
810
811
	list = read_reflections(infile);
	if ( list == NULL ) {
		ERROR("Couldn't read file '%s'\n", infile);
Thomas White's avatar
Thomas White committed
812
813
		return 1;
	}
Thomas White's avatar
Thomas White committed
814
	if ( check_list_symmetry(list, sym) ) {
Thomas White's avatar
Thomas White committed
815
		ERROR("The input reflection list does not appear to"
Thomas White's avatar
Thomas White committed
816
		      " have symmetry %s\n", symmetry_name(sym));
Thomas White's avatar
Thomas White committed
817
818
		return 1;
	}
Thomas White's avatar
Thomas White committed
819

820
	render_za(cell, list, boost, sym, wght, colscale,
821
	          rh, rk, rl, dh, dk, dl, outfile, scale_top, zone);
Thomas White's avatar
Thomas White committed
822

823
	free(pdb);
Thomas White's avatar
Thomas White committed
824
	free_symoplist(sym);
Thomas White's avatar
Thomas White committed
825
	reflist_free(list);
826
	if ( outfile != NULL ) free(outfile);
827

828
	return r;
829
}