From e9906650fd078b3a8c60581847d8d423f6c22e5e Mon Sep 17 00:00:00 2001
From: Thomas White <taw@physics.org>
Date: Thu, 15 Jun 2023 17:37:03 +0200
Subject: z-rotation gradients

---
 libcrystfel/src/crystfel-mille.c | 40 +++++++++++++++++++++++++++++++---------
 libcrystfel/src/detgeom.h        |  6 ++++++
 libcrystfel/src/predict-refine.c |  5 +++--
 libcrystfel/src/predict-refine.h |  3 +++
 4 files changed, 43 insertions(+), 11 deletions(-)

(limited to 'libcrystfel/src')

diff --git a/libcrystfel/src/crystfel-mille.c b/libcrystfel/src/crystfel-mille.c
index 88e6e981..6bdcb51a 100644
--- a/libcrystfel/src/crystfel-mille.c
+++ b/libcrystfel/src/crystfel-mille.c
@@ -62,9 +62,12 @@ int mille_label(int hierarchy_level, int member_index, char param)
 
 	label = 100000*hierarchy_level + 100*member_index;
 	switch ( param ) {
-		case 'x' : return label+1;
-		case 'y' : return label+2;
-		case 'z' : return label+3;
+		case 'x' : return label+1;  /* x-shift */
+		case 'y' : return label+2;  /* y-shift */
+		case 'z' : return label+3;  /* z-shift */
+		case 'a' : return label+4;  /* Rotation around x */
+		case 'b' : return label+5;  /* Rotation around y */
+		case 'c' : return label+6;  /* Rotation around z */
 		default : abort();
 	}
 }
@@ -78,12 +81,29 @@ void write_mille(Mille *mille, int n, UnitCell *cell,
 	float global_gradients[64];
 	int labels[6];
 	int i;
+	double asx, asy, asz, bsx, bsy, bsz, csx, csy, csz;
+
+	cell_get_reciprocal(cell, &asx, &asy, &asz,
+	                          &bsx, &bsy, &bsz,
+	                          &csx, &csy, &csz);
 
 	for ( i=0; i<n; i++ ) {
 
+		signed int h, k, l;
+		double xl, yl, zl, kpred;
+		double xpk, ypk;
 		int j;
 		const struct detgeom_panel_group *group;
 
+		get_indices(rps[i].refl, &h, &k, &l);
+		kpred = get_kpred(rps[i].refl);
+		xl = h*asx + k*bsx + l*csx;
+		yl = h*asy + k*bsy + l*csy;
+		zl = h*asz + k*bsz + l*csz;
+
+		twod_mapping(rps[i].peak->fs, rps[i].peak->ss, &xpk, &ypk,
+		             rps[i].panel, dx, dy);
+
 		/* x terms: local */
 		for ( j=0; j<9; j++ ) {
 			local_gradients[j] = x_gradient(rv[j], rps[i].refl,
@@ -99,12 +119,13 @@ void write_mille(Mille *mille, int n, UnitCell *cell,
 			labels[j] = mille_label(group->hierarchy_level, group->member_index, 'x');
 			j++;
 
-			global_gradients[j] = x_gradient(GPARAM_CLEN, rps[i].refl,
-			                                 cell, rps[i].panel);
+			global_gradients[j] = -xl / (kpred+zl);
 			labels[j] = mille_label(group->hierarchy_level, group->member_index, 'z');
 			j++;
 
-			/* FIXME: Rotations */
+			global_gradients[j] = -(ypk - group->cy);
+			labels[j] = mille_label(group->hierarchy_level, group->member_index, 'c');
+			j++;
 
 			group = group->parent;
 		}
@@ -130,12 +151,13 @@ void write_mille(Mille *mille, int n, UnitCell *cell,
 			labels[j] = mille_label(group->hierarchy_level, group->member_index, 'y');
 			j++;
 
-			global_gradients[j] = y_gradient(GPARAM_CLEN, rps[i].refl,
-			                                 cell, rps[i].panel);
+			global_gradients[j] = -yl / (kpred+zl);
 			labels[j] = mille_label(group->hierarchy_level, group->member_index, 'z');
 			j++;
 
-			/* FIXME: Rotations */
+			global_gradients[j] = xpk - group->cx;
+			labels[j] = mille_label(group->hierarchy_level, group->member_index, 'c');
+			j++;
 
 			group = group->parent;
 		}
diff --git a/libcrystfel/src/detgeom.h b/libcrystfel/src/detgeom.h
index 917b6bd3..22f5bd1a 100644
--- a/libcrystfel/src/detgeom.h
+++ b/libcrystfel/src/detgeom.h
@@ -98,6 +98,12 @@ struct detgeom_panel_group
 	int hierarchy_level;
 	int member_index;
 
+	/* Center of panel group, in lab coordinate system (metres)
+	 * This will be the rotation center. */
+	double cx;
+	double cy;
+	double cz;
+
 	/* If n_children > 0, here are the child groups */
 	struct detgeom_panel_group **children;
 
diff --git a/libcrystfel/src/predict-refine.c b/libcrystfel/src/predict-refine.c
index 4a53e8e0..5210846a 100644
--- a/libcrystfel/src/predict-refine.c
+++ b/libcrystfel/src/predict-refine.c
@@ -67,8 +67,9 @@ static const enum gparam rv[] =
 	GPARAM_DETY,
 };
 
-static void twod_mapping(double fs, double ss, double *px, double *py,
-                         struct detgeom_panel *p, double dx, double dy)
+
+void twod_mapping(double fs, double ss, double *px, double *py,
+                  struct detgeom_panel *p, double dx, double dy)
 {
 	double xs, ys;
 
diff --git a/libcrystfel/src/predict-refine.h b/libcrystfel/src/predict-refine.h
index 5e2552e6..d1e5fbe7 100644
--- a/libcrystfel/src/predict-refine.h
+++ b/libcrystfel/src/predict-refine.h
@@ -46,6 +46,9 @@ struct reflpeak {
  * Prediction refinement: refinement of indexing solutions before integration.
  */
 
+extern void twod_mapping(double fs, double ss, double *px, double *py,
+                         struct detgeom_panel *p, double dx, double dy);
+
 extern int refine_prediction(struct image *image, Crystal *cr, Mille *mille);
 
 extern int refine_radius(Crystal *cr, struct image *image);
-- 
cgit v1.2.3