process_hkl: Merge with symmetry

src/process_hkl.c

@@ -26,6 +26,7 @@
 #include "sfac.h"
 #include "reflections.h"
 #include "likelihood.h"
+#include "symmetry.h"
 
 
 /* Number of divisions for R vs |q| graphs */
@@ -188,21 +189,25 @@ static void process_reflections(double *ref, unsigned int *counts,
 
 static void merge_pattern(double *model, const double *new,
                           unsigned int *model_counts,
-                          ReflItemList *items, int mo, int sum)
+                          ReflItemList *items, int mo, int sum,
+                          const char *symm)
 {
 	int i;
 
 	for ( i=0; i<num_items(items); i++ ) {
 
 		double intensity;
+		signed int hs, ks, ls;
 		signed int h, k, l;
 		struct refl_item *item;
 
 		item = get_item(items, i);
 
-		h = item->h;
-		k = item->k;
-		l = item->l;
+		hs = item->h;
+		ks = item->k;
+		ls = item->l;
+
+		get_asymm(hs, ks, ls, &h, &k, &l, symm);
 
 		intensity = lookup_intensity(new, h, k, l);
 
@@ -268,7 +273,7 @@ int main(int argc, char *argv[])
 		{"sum",                0, &config_sum,         1},
 		{"detwin",             0, &config_detwin,      1},
 		{"scale",              0, &config_scale,       1},
-		{"symmetry",           0, NULL,               'y'},
+		{"symmetry",           1, NULL,               'y'},
 		{"pdb",                1, NULL,               'p'},
 		{0, 0, NULL, 0}
 	};
@@ -413,7 +418,7 @@ int main(int argc, char *argv[])
 
 			/* Start of second or later pattern */
 			merge_pattern(model, new_pattern, model_counts,
-			              items, config_maxonly, config_sum);
+			              items, config_maxonly, config_sum, sym);
 
 			if ( (trueref != NULL) && config_every
 			    && (n_patterns % config_every == 0) ) {

src/symmetry.c

@@ -18,6 +18,79 @@
 #include <stdio.h>
 #include <math.h>
 
-void apply_symmetry(double ref, unsigned int *counts, const char *sym)
+#include "utils.h"
+
+
+/* Conditions for a reflection to be in the asymmetric unit cell */
+#define COND_1(h, k, l) (1)
+#define COND_6MMM(h, k, i, l) ( (h>=0) && (k>=0) && (l>=0) && ((h>k)||((h==0)&&(k==0))) )
+
+
+/* Macros for checking the above conditions and returning if satisfied */
+#define CHECK_COND_FOURIDX(h, k, i, l, cond)       \
+	if ( COND_##cond((h), (k), (i), (l)) ) {   \
+		*hp = (h);  *kp = (k);  *lp = (l); \
+		return;                            \
+	}
+
+#define CHECK_COND_THREEIDX(h, k, l, cond)          \
+	if ( COND_##cond((h), (k), (l)) ) {         \
+		*hp = (h);  *kp = (k);  *lp = (l);  \
+		return;                             \
+	}
+
+
+/* Abort macro if no match found */
+#define SYM_ABORT                                             \
+	ERROR("No match in %s for %i %i %i\n", sym, h, k, l); \
+	abort();
+
+
+void get_asymm(signed int h, signed int k, signed int l,
+               signed int *hp, signed int *kp, signed int *lp,
+               const char *sym)
 {
+	if ( strcmp(sym, "1") == 0 ) {
+		CHECK_COND_THREEIDX(h, k, l, 1);
+		SYM_ABORT;
+	}
+
+	if ( strcmp(sym, "6/mmm") == 0 ) {
+
+		const signed int i = h+k;
+
+		CHECK_COND_FOURIDX(h, k, i, l, 6MMM);
+		CHECK_COND_FOURIDX(h, i, k, l, 6MMM);
+		CHECK_COND_FOURIDX(k, h, i, l, 6MMM);
+		CHECK_COND_FOURIDX(k, i, h, l, 6MMM);
+		CHECK_COND_FOURIDX(i, h, k, l, 6MMM);
+		CHECK_COND_FOURIDX(i, k, h, l, 6MMM);
+
+		CHECK_COND_FOURIDX(h, k, i, -l, 6MMM);
+		CHECK_COND_FOURIDX(h, i, k, -l, 6MMM);
+		CHECK_COND_FOURIDX(k, h, i, -l, 6MMM);
+		CHECK_COND_FOURIDX(k, i, h, -l, 6MMM);
+		CHECK_COND_FOURIDX(i, h, k, -l, 6MMM);
+		CHECK_COND_FOURIDX(i, k, h, -l, 6MMM);
+
+		CHECK_COND_FOURIDX(-h, -k, -i, l, 6MMM);
+		CHECK_COND_FOURIDX(-h, -i, -k, l, 6MMM);
+		CHECK_COND_FOURIDX(-k, -h, -i, l, 6MMM);
+		CHECK_COND_FOURIDX(-k, -i, -h, l, 6MMM);
+		CHECK_COND_FOURIDX(-i, -h, -k, l, 6MMM);
+		CHECK_COND_FOURIDX(-i, -k, -h, l, 6MMM);
+
+		CHECK_COND_FOURIDX(-h, -k, -i, -l, 6MMM);
+		CHECK_COND_FOURIDX(-h, -i, -k, -l, 6MMM);
+		CHECK_COND_FOURIDX(-k, -h, -i, -l, 6MMM);
+		CHECK_COND_FOURIDX(-k, -i, -h, -l, 6MMM);
+		CHECK_COND_FOURIDX(-i, -h, -k, -l, 6MMM);
+		CHECK_COND_FOURIDX(-i, -k, -h, -l, 6MMM);
+
+		SYM_ABORT;  /* Should never reach here */
+
+	}
+
+	ERROR("Unknown point group '%s'\n", sym);
+	abort();
 }

src/symmetry.h

@@ -17,6 +17,9 @@
 #include <config.h>
 #endif
 
-extern void apply_symmetry(double ref, unsigned int *counts, const char *sym);
+
+extern void get_asymm(signed int h, signed int k, signed int l,
+                      signed int *hp, signed int *kp, signed int *lp,
+                      const char *sym);
 
 #endif	/* SYMMETRY_H */