diff options
-rw-r--r-- | libcrystfel/src/geometry.c | 19 | ||||
-rw-r--r-- | libcrystfel/src/geometry.h | 2 | ||||
-rw-r--r-- | src/partialator.c | 2 |
3 files changed, 17 insertions, 6 deletions
diff --git a/libcrystfel/src/geometry.c b/libcrystfel/src/geometry.c index 3586ed0b..21f81915 100644 --- a/libcrystfel/src/geometry.c +++ b/libcrystfel/src/geometry.c @@ -105,7 +105,7 @@ static signed int locate_peak(double x, double y, double z, double k, static double partiality(PartialityModel pmodel, double rlow, double rmid, double rhigh, - double r) + double r, double D) { double qlow, qhigh; double plow, phigh; @@ -135,6 +135,11 @@ static double partiality(PartialityModel pmodel, case PMODEL_THIN: return 1.0 - (rmid*rmid)/(r*r); + case PMODEL_SCSPHERE: + plow = 3.0*qlow*qlow - 2.0*qlow*qlow*qlow; + phigh = 3.0*qhigh*qhigh - 2.0*qhigh*qhigh*qhigh; + return 4.0*(plow-phigh)*r / (3.0*D); + } } @@ -153,7 +158,7 @@ static Reflection *check_reflection(struct image *image, Crystal *cryst, Reflection *refl; double cet, cez; /* Centre of Ewald sphere */ double pr; - double L; + double L, D; double del; /* Don't predict 000 */ @@ -197,10 +202,10 @@ static Reflection *check_reflection(struct image *image, Crystal *cryst, /* Conditions for reflection to be excited at all */ switch ( pmodel ) { - default: case PMODEL_UNITY: /* PMODEL_UNITY shouldn't end up here */ case PMODEL_SPHERE: case PMODEL_GAUSSIAN: + case PMODEL_SCSPHERE: if ( (signbit(rlow) == signbit(rhigh)) && (fabs(rlow) > pr) && (fabs(rhigh) > pr) ) return NULL; @@ -212,19 +217,21 @@ static Reflection *check_reflection(struct image *image, Crystal *cryst, } + D = rlow - rhigh; + /* Lorentz factor is determined direction from the r values, before * clamping. The multiplication by 0.01e9 to make the * correction factor vaguely near 1. */ switch ( pmodel ) { - default: case PMODEL_SPHERE: case PMODEL_GAUSSIAN: - L = LORENTZ_SCALE / (rlow - rhigh); + L = LORENTZ_SCALE / D; break; case PMODEL_UNITY: /* PMODEL_UNITY shouldn't end up here */ case PMODEL_THIN: + case PMODEL_SCSPHERE: L = 1.0; break; @@ -257,7 +264,7 @@ static Reflection *check_reflection(struct image *image, Crystal *cryst, } /* Calculate partiality */ - part = partiality(pmodel, rlow, rmid, rhigh, pr); + part = partiality(pmodel, rlow, rmid, rhigh, pr, D); /* Add peak to list */ refl = reflection_new(h, k, l); diff --git a/libcrystfel/src/geometry.h b/libcrystfel/src/geometry.h index 8041936a..d8d226f0 100644 --- a/libcrystfel/src/geometry.h +++ b/libcrystfel/src/geometry.h @@ -51,6 +51,7 @@ extern "C" { * @PMODEL_UNITY : Set all all partialities and Lorentz factors to 1. * @PMODEL_GAUSSIAN : Gaussian profiles in 3D * @PMODEL_THIN : Thin Ewald sphere intersecting sphere + * @PMODEL_SCSPHERE : Sphere model with source coverage factor included * * A %PartialityModel describes a geometrical model which can be used to * calculate spot partialities and Lorentz correction factors. @@ -61,6 +62,7 @@ typedef enum { PMODEL_UNITY, PMODEL_GAUSSIAN, PMODEL_THIN, + PMODEL_SCSPHERE, } PartialityModel; diff --git a/src/partialator.c b/src/partialator.c index e5083611..0d532997 100644 --- a/src/partialator.c +++ b/src/partialator.c @@ -350,6 +350,8 @@ int main(int argc, char *argv[]) pmodel = PMODEL_GAUSSIAN; } else if ( strcmp(pmodel_str, "thin") == 0 ) { pmodel = PMODEL_THIN; + } else if ( strcmp(pmodel_str, "scsphere") == 0 ) { + pmodel = PMODEL_SCSPHERE; } else { ERROR("Unknown partiality model '%s'.\n", pmodel_str); return 1; |