diff options
| -rw-r--r-- | libcrystfel/src/geometry.c | 65 | ||||
| -rw-r--r-- | libcrystfel/src/geometry.h | 1 | ||||
| -rw-r--r-- | src/partialator.c | 2 |
3 files changed, 68 insertions, 0 deletions
diff --git a/libcrystfel/src/geometry.c b/libcrystfel/src/geometry.c index 4f9515c5..e5f33d64 100644 --- a/libcrystfel/src/geometry.c +++ b/libcrystfel/src/geometry.c @@ -711,6 +711,67 @@ static void ginn_spectrum_partialities(Crystal *cryst) } +static void ewald_offset_partialities(Crystal *cryst) +{ + RefList *list; + Reflection *refl; + RefListIterator *iter; + double r0, m; + struct image *image; + UnitCell *cell; + double asx, asy, asz, bsx, bsy, bsz, csx, csy, csz; + + list = crystal_get_reflections(cryst); + if ( list == NULL ) { + ERROR("No reflections for partiality calculation!\n"); + return; + } + + image = crystal_get_image(cryst); + if ( image == NULL ) { + ERROR("No image for partiality calculation!\n"); + return; + } + + cell = crystal_get_cell(cryst); + if ( cell == NULL ) { + ERROR("No unit cell for partiality calculation!\n"); + return; + } + cell_get_reciprocal(cell, &asx, &asy, &asz, + &bsx, &bsy, &bsz, + &csx, &csy, &csz); + + r0 = fabs(crystal_get_profile_radius(cryst)); + m = crystal_get_mosaicity(cryst); + + for ( refl = first_refl(crystal_get_reflections(cryst), &iter); + refl != NULL; + refl = next_refl(refl, iter) ) + { + signed int h, k, l; + double xl, yl, zl; + double q2, R, t; + + get_symmetric_indices(refl, &h, &k, &l); + xl = h*asx + k*bsx + l*csx; + yl = h*asy + k*bsy + l*csy; + zl = h*asz + k*bsz + l*csz; + + /* Radius of rlp profile */ + q2 = xl*xl + yl*yl + zl*zl; + R = r0 + m * sqrt(q2); + + /* Excitation error */ + t = get_exerr(refl); + + set_partiality(refl, exp(-(t*t)/(R*R))); + set_lorentz(refl, 1.0); + + } +} + + /** * \param cryst A \ref Crystal * \param pmodel A \ref PartialityModel @@ -735,6 +796,10 @@ void calculate_partialities(Crystal *cryst, PartialityModel pmodel) ginn_spectrum_partialities(cryst); break; + case PMODEL_OFFSET : + ewald_offset_partialities(cryst); + break; + case PMODEL_RANDOM : set_random_partialities(cryst); break; diff --git a/libcrystfel/src/geometry.h b/libcrystfel/src/geometry.h index a9000428..e6e862a5 100644 --- a/libcrystfel/src/geometry.h +++ b/libcrystfel/src/geometry.h @@ -59,6 +59,7 @@ typedef enum { PMODEL_UNITY, /**< Set all partialities and Lorentz factors to 1. */ PMODEL_XSPHERE, /**< Flat sphere model with super-Gaussian spectrum */ + PMODEL_OFFSET, /**< Ewald offset model for monochromatic beam */ PMODEL_RANDOM, /**< Randomly assigned partialities */ } PartialityModel; diff --git a/src/partialator.c b/src/partialator.c index 5543b041..68c886c9 100644 --- a/src/partialator.c +++ b/src/partialator.c @@ -1267,6 +1267,8 @@ int main(int argc, char *argv[]) pmodel = PMODEL_UNITY; } else if ( strcmp(pmodel_str, "xsphere") == 0 ) { pmodel = PMODEL_XSPHERE; + } else if ( strcmp(pmodel_str, "offset") == 0 ) { + pmodel = PMODEL_OFFSET; } else if ( strcmp(pmodel_str, "random") == 0 ) { pmodel = PMODEL_RANDOM; } else { |