aboutsummaryrefslogtreecommitdiff
path: root/libcrystfel/src
diff options
context:
space:
mode:
Diffstat (limited to 'libcrystfel/src')
-rw-r--r--libcrystfel/src/cell-utils.c97
-rw-r--r--libcrystfel/src/cell-utils.h7
-rw-r--r--libcrystfel/src/crystfel-mille.c30
-rw-r--r--libcrystfel/src/datatemplate.c61
-rw-r--r--libcrystfel/src/detgeom.c1
-rw-r--r--libcrystfel/src/geometry.c3
-rw-r--r--libcrystfel/src/image-hdf5.c77
-rw-r--r--libcrystfel/src/image.c3
-rw-r--r--libcrystfel/src/index.c7
-rw-r--r--libcrystfel/src/indexers/asdf.c2
-rw-r--r--libcrystfel/src/indexers/pinkindexer.c1
-rw-r--r--libcrystfel/src/indexers/taketwo.c82
-rw-r--r--libcrystfel/src/integration.c7
-rw-r--r--libcrystfel/src/peakfinder8.c1
-rw-r--r--libcrystfel/src/predict-refine.c31
-rw-r--r--libcrystfel/src/utils.c26
-rw-r--r--libcrystfel/src/utils.h1
17 files changed, 301 insertions, 136 deletions
diff --git a/libcrystfel/src/cell-utils.c b/libcrystfel/src/cell-utils.c
index 80cdd955..6e5d65b7 100644
--- a/libcrystfel/src/cell-utils.c
+++ b/libcrystfel/src/cell-utils.c
@@ -3,12 +3,12 @@
*
* Unit Cell utility functions
*
- * Copyright © 2012-2021 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2012-2024 Deutsches Elektronen-Synchrotron DESY,
* a research centre of the Helmholtz Association.
* Copyright © 2012 Lorenzo Galli
*
* Authors:
- * 2009-2021 Thomas White <taw@physics.org>
+ * 2009-2024 Thomas White <taw@physics.org>
* 2012 Lorenzo Galli
*
* This file is part of CrystFEL.
@@ -42,6 +42,7 @@
#include "cell-utils.h"
#include "utils.h"
#include "image.h"
+#include "symmetry.h"
/**
@@ -2245,9 +2246,9 @@ static IntegerMatrix *check_permutations(UnitCell *cell_reduced, UnitCell *refer
best_diff = diff;
n_best = 0;
- }
+ best_m[n_best++] = intmat_copy(m);
- if ( diff < 1.001*best_diff ) {
+ } else if ( diff < 1.001*best_diff ) {
/* If the new solution is the same as the
* previous one, add it to the list */
@@ -2414,3 +2415,91 @@ UnitCell *compare_reindexed_cell_parameters(UnitCell *cell_in,
return match;
}
+
+
+static char *add_unique_axis(const char *inp, char ua)
+{
+ char *pg = cfmalloc(64);
+ if ( pg == NULL ) return NULL;
+ snprintf(pg, 63, "%s_ua%c", inp, ua);
+ return pg;
+}
+
+
+static char *get_chiral_holohedry(UnitCell *cell)
+{
+ LatticeType lattice = cell_get_lattice_type(cell);
+ char *pg;
+ int add_ua = 1;
+
+ switch (lattice)
+ {
+ case L_TRICLINIC:
+ pg = "1";
+ add_ua = 0;
+ break;
+
+ case L_MONOCLINIC:
+ pg = "2";
+ break;
+
+ case L_ORTHORHOMBIC:
+ pg = "222";
+ add_ua = 0;
+ break;
+
+ case L_TETRAGONAL:
+ pg = "422";
+ break;
+
+ case L_RHOMBOHEDRAL:
+ pg = "3_R";
+ add_ua = 0;
+ break;
+
+ case L_HEXAGONAL:
+ if ( cell_get_centering(cell) == 'H' ) {
+ pg = "3_H";
+ add_ua = 0;
+ } else {
+ pg = "622";
+ }
+ break;
+
+ case L_CUBIC:
+ pg = "432";
+ add_ua = 0;
+ break;
+
+ default:
+ pg = "error";
+ break;
+ }
+
+ if ( add_ua ) {
+ return add_unique_axis(pg, cell_get_unique_axis(cell));
+ } else {
+ return cfstrdup(pg);
+ }
+}
+
+
+/**
+ * \param cell: A UnitCell
+ *
+ * Determines the chiral point group of the lattice, based on the lattice type
+ * and unique axis of \p cell (the lattice parameters are not checked).
+ *
+ * This point group is often needed within indexing algorithms.
+ *
+ * \returns a newly allocated SymOpList.
+ *
+ */
+SymOpList *get_lattice_symmetry(UnitCell *cell)
+{
+ SymOpList *rawList;
+ char *pg = get_chiral_holohedry(cell);
+ rawList = get_pointgroup(pg);
+ cffree(pg);
+ return rawList;
+}
diff --git a/libcrystfel/src/cell-utils.h b/libcrystfel/src/cell-utils.h
index b473a4c5..de200bcb 100644
--- a/libcrystfel/src/cell-utils.h
+++ b/libcrystfel/src/cell-utils.h
@@ -3,12 +3,12 @@
*
* Unit Cell utility functions
*
- * Copyright © 2012-2021 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2012-2024 Deutsches Elektronen-Synchrotron DESY,
* a research centre of the Helmholtz Association.
* Copyright © 2012 Lorenzo Galli
*
* Authors:
- * 2009-2020 Thomas White <taw@physics.org>
+ * 2009-2024 Thomas White <taw@physics.org>
* 2012 Lorenzo Galli
*
* This file is part of CrystFEL.
@@ -34,6 +34,7 @@
#include <gsl/gsl_matrix.h>
#include "cell.h"
+#include "symmetry.h"
#ifdef __cplusplus
extern "C" {
@@ -102,6 +103,8 @@ extern UnitCell *compare_reindexed_cell_parameters(UnitCell *cell_in,
const double *tols,
RationalMatrix **pmb);
+extern SymOpList *get_lattice_symmetry(UnitCell *cell);
+
#ifdef __cplusplus
}
#endif
diff --git a/libcrystfel/src/crystfel-mille.c b/libcrystfel/src/crystfel-mille.c
index 72f89425..c2fe2a57 100644
--- a/libcrystfel/src/crystfel-mille.c
+++ b/libcrystfel/src/crystfel-mille.c
@@ -72,6 +72,10 @@ struct mille
int *int_arr;
int max_entries;
int n;
+
+ int *have_local;
+ int n_local;
+
FILE *fh;
};
@@ -115,6 +119,17 @@ static void mille_add_measurement(Mille *m,
m->max_entries = new_max_entries;
}
+ if ( NLC > m->n_local ) {
+ int i;
+ int *new_have_local = cfrealloc(m->have_local, NLC*sizeof(int));
+ if ( new_have_local == NULL ) return;
+ m->have_local = new_have_local;
+ for ( i=m->n_local; i<NLC; i++ ) {
+ m->have_local[i] = 0;
+ }
+ m->n_local = NLC;
+ }
+
/* The measurement */
m->float_arr[m->n] = rMeas;
m->int_arr[m->n] = 0;
@@ -126,6 +141,7 @@ static void mille_add_measurement(Mille *m,
m->float_arr[m->n] = derLc[i];
m->int_arr[m->n] = i+1;
m->n++;
+ m->have_local[i] = 1;
}
}
@@ -292,6 +308,8 @@ Mille *crystfel_mille_new(const char *outFileName)
m->n = 0;
m->float_arr = NULL;
m->int_arr = NULL;
+ m->have_local = NULL;
+ m->n_local = 0;
m->fh = fopen(outFileName, "wb");
if ( m->fh == NULL ) {
@@ -318,11 +336,15 @@ void crystfel_mille_free(Mille *m)
void crystfel_mille_delete_last_record(Mille *m)
{
m->n = 0;
+ cffree(m->have_local);
+ m->have_local = NULL;
+ m->n_local = 0;
}
void crystfel_mille_write_record(Mille *m)
{
+ int i;
float nf = 0.0;
int ni = 0;
int nw = (m->n * 2)+2;
@@ -330,6 +352,14 @@ void crystfel_mille_write_record(Mille *m)
/* Don't write empty records */
if ( m->n == 0 ) return;
+ /* Don't write records with incomplete local gradients */
+ for ( i=0; i<m->n_local; i++ ) {
+ if ( !m->have_local[i] ) {
+ crystfel_mille_delete_last_record(m);
+ return;
+ }
+ }
+
fwrite(&nw, sizeof(int), 1, m->fh);
fwrite(&nf, sizeof(float), 1, m->fh);
diff --git a/libcrystfel/src/datatemplate.c b/libcrystfel/src/datatemplate.c
index a7b8950f..c97bd6d7 100644
--- a/libcrystfel/src/datatemplate.c
+++ b/libcrystfel/src/datatemplate.c
@@ -3,11 +3,11 @@
*
* Data template structure
*
- * Copyright © 2019-2021 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2019-2024 Deutsches Elektronen-Synchrotron DESY,
* a research centre of the Helmholtz Association.
*
* Authors:
- * 2019-2021 Thomas White <taw@physics.org>
+ * 2019-2024 Thomas White <taw@physics.org>
*
* This file is part of CrystFEL.
*
@@ -830,12 +830,35 @@ static int parse_peak_layout(const char *val,
return 1;
}
+#define MAX_FOR_LATER (1024)
+
+struct forlater
+{
+ char *keys[MAX_FOR_LATER];
+ char *vals[MAX_FOR_LATER];
+ int n_forlater;
+};
+
+
+static void store_for_later(struct forlater *fl, const char *key, const char *val)
+{
+ if ( fl->n_forlater >= MAX_FOR_LATER ) {
+ ERROR("Too many lines stored.\n");
+ return;
+ }
+
+ fl->keys[fl->n_forlater] = strdup(key);
+ fl->vals[fl->n_forlater] = strdup(val);
+ fl->n_forlater++;
+}
+
static int parse_toplevel(DataTemplate *dt,
const char *key,
const char *val,
struct panel_template *defaults,
- int *defaults_updated)
+ int *defaults_updated,
+ struct forlater *for_later)
{
if ( strcmp(key, "detector_shift_x") == 0 ) {
dt->shift_x_from = cfstrdup(val);
@@ -883,8 +906,12 @@ static int parse_toplevel(DataTemplate *dt,
} else if ( strncmp(key, "group_", 6) == 0 ) {
- if ( parse_group(key+6, dt, val) ) {
- return 1;
+ if ( for_later != NULL ) {
+ store_for_later(for_later, key, val);
+ } else {
+ if ( parse_group(key+6, dt, val) ) {
+ return 1;
+ }
}
} else {
@@ -1042,6 +1069,9 @@ DataTemplate *data_template_new_from_string(const char *string_in)
size_t len;
struct panel_template defaults;
int have_unused_defaults = 0;
+ struct forlater lines_for_later;
+
+ lines_for_later.n_forlater = 0;
dt = cfcalloc(1, sizeof(DataTemplate));
if ( dt == NULL ) return NULL;
@@ -1179,10 +1209,10 @@ DataTemplate *data_template_new_from_string(const char *string_in)
/* Top-level option */
if ( parse_toplevel(dt, line, val,
&defaults,
- &have_unused_defaults) )
+ &have_unused_defaults,
+ &lines_for_later) )
{
- ERROR("Invalid top-level line '%s'\n",
- line);
+ ERROR("Invalid top-level line '%s'\n", line);
reject = 1;
}
cffree(line);
@@ -1217,6 +1247,21 @@ DataTemplate *data_template_new_from_string(const char *string_in)
} while ( !done );
+ for ( i=0; i<lines_for_later.n_forlater; i++ ) {
+ char *key = lines_for_later.keys[i];
+ char *val = lines_for_later.vals[i];
+ if ( parse_toplevel(dt, key, val,
+ &defaults,
+ &have_unused_defaults,
+ NULL) )
+ {
+ ERROR("Invalid top-level line '%s' = '%s'\n", key, val);
+ reject = 1;
+ }
+ cffree(key);
+ cffree(val);
+ }
+
if ( dt->n_panels == 0 ) {
ERROR("No panel descriptions in geometry file.\n");
cffree(dt);
diff --git a/libcrystfel/src/detgeom.c b/libcrystfel/src/detgeom.c
index 39b00663..2b38664c 100644
--- a/libcrystfel/src/detgeom.c
+++ b/libcrystfel/src/detgeom.c
@@ -269,6 +269,7 @@ gsl_matrix **make_panel_minvs(struct detgeom *dg)
gsl_matrix_set(M, 2, 2, p->pixel_pitch*p->ssz);
Minvs[i] = matrix_invert(M);
+ gsl_matrix_free(M);
if ( Minvs[i] == NULL ) {
ERROR("Failed to calculate inverse panel matrix for %s\n",
p->name);
diff --git a/libcrystfel/src/geometry.c b/libcrystfel/src/geometry.c
index 10a7fa19..7bcb8d8f 100644
--- a/libcrystfel/src/geometry.c
+++ b/libcrystfel/src/geometry.c
@@ -306,7 +306,7 @@ static Reflection *check_reflection(struct image *image, Crystal *cryst,
z -= g.kcen;
/* Three because the general case fails in extreme cases */
- if ( w0 / w1 <= DBL_MIN ) {
+ if ( isinf(w0) || (w0 / w1 <= DBL_MIN) ) {
/* 'Laue' corner case */
kpred = g.kcen;
@@ -352,6 +352,7 @@ static Reflection *check_reflection(struct image *image, Crystal *cryst,
/* Revert the 'Lorentz' factor */
partiality *= sqrt( ( R*R + M2_k/sumw_k) / ( R*R ) );
+ if ( isnan(partiality) ) partiality = 0.0;
if ( (updateme == NULL) && ( partiality < min_partiality ) ) return NULL;
diff --git a/libcrystfel/src/image-hdf5.c b/libcrystfel/src/image-hdf5.c
index a148c307..93ea2d8d 100644
--- a/libcrystfel/src/image-hdf5.c
+++ b/libcrystfel/src/image-hdf5.c
@@ -414,6 +414,7 @@ static int load_hdf5_hyperslab(struct panel_template *p,
if ( ndims < 0 ) {
ERROR("Failed to get number of dimensions for panel %s\n",
p->name);
+ H5Sclose(dataspace);
H5Dclose(dh);
return 1;
}
@@ -436,6 +437,7 @@ static int load_hdf5_hyperslab(struct panel_template *p,
p->name, ndims, total_dt_dims,
total_dt_dims - plh_dt_dims);
H5Dclose(dh);
+ H5Sclose(dataspace);
return 1;
}
} else {
@@ -453,6 +455,7 @@ static int load_hdf5_hyperslab(struct panel_template *p,
cffree(f_offset);
cffree(f_count);
H5Dclose(dh);
+ H5Sclose(dataspace);
return 1;
}
@@ -502,6 +505,7 @@ static int load_hdf5_hyperslab(struct panel_template *p,
cffree(f_offset);
cffree(f_count);
H5Dclose(dh);
+ H5Sclose(dataspace);
return 1;
}
@@ -511,6 +515,8 @@ static int load_hdf5_hyperslab(struct panel_template *p,
profile_start("H5Dread");
r = H5Dread(dh, el_type, memspace, dataspace, H5P_DEFAULT, data);
+ H5Sclose(memspace);
+ H5Sclose(dataspace);
profile_end("H5Dread");
if ( r < 0 ) {
ERROR("Couldn't read data for panel %s\n",
@@ -538,6 +544,7 @@ static int load_hdf5_hyperslab(struct panel_template *p,
static hid_t open_hdf5_file(const char *filename)
{
hid_t fh;
+ hid_t fapl;
if ( access(filename, R_OK) == -1 ) {
ERROR("File does not exist or cannot be read: %s\n",
@@ -545,7 +552,10 @@ static hid_t open_hdf5_file(const char *filename)
return -1;
}
- fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT);
+ fapl = H5Pcreate(H5P_FILE_ACCESS);
+ H5Pset_fclose_degree(fapl, H5F_CLOSE_STRONG);
+ fh = H5Fopen(filename, H5F_ACC_RDONLY, fapl);
+ H5Pclose(fapl);
if ( fh < 0 ) {
ERROR("Couldn't open HDF5 file: %s\n", filename);
return -1;
@@ -650,6 +660,7 @@ int image_hdf5_read_satmap(struct panel_template *p,
sizeof(float), 1, map_location, NULL) )
{
ERROR("Failed to load saturation map data\n");
+ close_hdf5(fh);
return 1;
}
@@ -687,6 +698,7 @@ int image_hdf5_read_mask(struct panel_template *p,
sizeof(int), 1, mask_location, NULL) )
{
ERROR("Failed to load mask data\n");
+ close_hdf5(fh);
cffree(mask);
return 1;
}
@@ -718,6 +730,7 @@ static char *read_single_fixed_string(hid_t dh)
sh = H5Dget_space(dh);
if ( H5Sget_simple_extent_ndims(sh) ) {
ERROR("Non-scalar string\n");
+ H5Sclose(sh);
return NULL;
}
@@ -727,6 +740,7 @@ static char *read_single_fixed_string(hid_t dh)
tmp = cfmalloc(size+1);
if ( tmp == NULL ) {
H5Tclose(type);
+ H5Sclose(sh);
return NULL;
}
r = H5Dread(dh, type, sh, H5S_ALL, H5P_DEFAULT, tmp);
@@ -832,6 +846,7 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
ndims = H5Sget_simple_extent_ndims(sh);
if ( ndims > 64 ) {
ERROR("Too many dimensions for numeric value\n");
+ H5Sclose(sh);
close_hdf5(fh);
cffree(subst_name);
return 1;
@@ -856,11 +871,15 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
ERROR("Couldn't read scalar value from %s.\n",
subst_name);
cffree(subst_name);
+ H5Sclose(sh);
+ H5Sclose(ms);
close_hdf5(fh);
return 1;
}
image_cache_header_float(image, name, val);
cffree(subst_name);
+ H5Sclose(sh);
+ H5Sclose(ms);
return 0;
} else if ( class == H5T_INTEGER ) {
@@ -868,6 +887,8 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
int val;
r = H5Dread(dh, H5T_NATIVE_INT, ms, sh, H5P_DEFAULT,
&val);
+ H5Sclose(sh);
+ H5Sclose(ms);
if ( r < 0 ) {
ERROR("Couldn't read scalar value from %s.\n",
subst_name);
@@ -912,12 +933,16 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
cffree(subst_name);
close_hdf5(fh);
+ H5Sclose(sh);
+ H5Sclose(ms);
return rv;
} else {
/* Should never be reached */
ERROR("Invalid HDF5 class %i\n", class);
cffree(subst_name);
+ H5Sclose(sh);
+ H5Sclose(ms);
return 1;
}
}
@@ -927,6 +952,8 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
ERROR("Couldn't parse event '%s'\n");
close_hdf5(fh);
cffree(subst_name);
+ H5Sclose(sh);
+ H5Sclose(ms);
return 1;
}
@@ -936,6 +963,8 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
ERROR("Couldn't allocate dimension arrays\n");
close_hdf5(fh);
cffree(subst_name);
+ H5Sclose(sh);
+ H5Sclose(ms);
return 1;
}
@@ -953,6 +982,8 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
subst_name, i,
dim_vals[dim_val_pos], size[i]);
close_hdf5(fh);
+ H5Sclose(sh);
+ H5Sclose(ms);
cffree(subst_name);
cffree(dim_vals);
return 1;
@@ -979,6 +1010,8 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
cffree(f_offset);
cffree(f_count);
cffree(subst_name);
+ H5Sclose(sh);
+ H5Sclose(ms);
close_hdf5(fh);
return 1;
}
@@ -986,12 +1019,13 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
cffree(f_offset);
cffree(f_count);
- ms = H5Screate_simple(1,msdims,NULL);
check = H5Sselect_hyperslab(ms, H5S_SELECT_SET,
m_offset, NULL, m_count, NULL);
if ( check < 0 ) {
ERROR("Error selecting memory dataspace for header value\n");
close_hdf5(fh);
+ H5Sclose(sh);
+ H5Sclose(ms);
cffree(subst_name);
return 1;
}
@@ -1000,6 +1034,8 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
double val;
r = H5Dread(dh, H5T_NATIVE_DOUBLE, ms, sh, H5P_DEFAULT, &val);
+ H5Sclose(sh);
+ H5Sclose(ms);
if ( r < 0 ) {
ERROR("Couldn't read value.\n");
close_hdf5(fh);
@@ -1016,6 +1052,8 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
int val;
r = H5Dread(dh, H5T_NATIVE_INT, ms, sh, H5P_DEFAULT, &val);
+ H5Sclose(sh);
+ H5Sclose(ms);
if ( r < 0 ) {
ERROR("Couldn't read value.\n");
close_hdf5(fh);
@@ -1041,6 +1079,8 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
char *val;
rv = H5Dread(dh, stype, ms, sh, H5P_DEFAULT, &val);
+ H5Sclose(sh);
+ H5Sclose(ms);
if ( rv < 0 ) {
ERROR("Can't read HDF5 vlen string from array - %s\n",
subst_name);
@@ -1069,10 +1109,14 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
val = cfmalloc(ssize+1);
if ( val == NULL ) {
close_hdf5(fh);
+ H5Sclose(ms);
+ H5Sclose(sh);
cffree(subst_name);
return 1;
}
rv = H5Dread(dh, stype, ms, sh, H5P_DEFAULT, val);
+ H5Sclose(sh);
+ H5Sclose(ms);
H5Tclose(stype);
if ( rv < 0 ) {
ERROR("Couldn't read HDF5 fixed string from array - %s\n",
@@ -1097,6 +1141,8 @@ int image_hdf5_read_header_to_cache(struct image *image, const char *name)
} else {
/* Should never be reached */
ERROR("Invalid HDF5 class %i\n", class);
+ H5Sclose(sh);
+ H5Sclose(ms);
close_hdf5(fh);
cffree(subst_name);
return 1;
@@ -1303,7 +1349,7 @@ ImageFeatureList *image_hdf5_read_peaks_cxi(const DataTemplate *dtempl,
int half_pixel_shift)
{
ImageFeatureList *features;
- hid_t fh;
+ hid_t fh, fapl;
char path_n[1024];
char path_x[1024];
char path_y[1024];
@@ -1355,7 +1401,10 @@ ImageFeatureList *image_hdf5_read_peaks_cxi(const DataTemplate *dtempl,
snprintf(path_y, 1024, "%s/peakYPosRaw", subst_name);
snprintf(path_i, 1024, "%s/peakTotalIntensity", subst_name);
- fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT);
+ fapl = H5Pcreate(H5P_FILE_ACCESS);
+ H5Pset_fclose_degree(fapl, H5F_CLOSE_STRONG);
+ fh = H5Fopen(filename, H5F_ACC_RDONLY, fapl);
+ H5Pclose(fapl);
if ( fh < 0 ) {
ERROR("Couldn't open file (peaks/cxi): %s\n", filename);
cffree(subst_name);
@@ -1431,7 +1480,7 @@ ImageFeatureList *image_hdf5_read_peaks_hdf5(const DataTemplate *dtempl,
const char *event,
int half_pixel_shift)
{
- hid_t fh, dh, sh;
+ hid_t fh, dh, sh, fapl;
hsize_t size[2];
hsize_t max_size[2];
int i;
@@ -1453,7 +1502,10 @@ ImageFeatureList *image_hdf5_read_peaks_hdf5(const DataTemplate *dtempl,
return NULL;
}
- fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT);
+ fapl = H5Pcreate(H5P_FILE_ACCESS);
+ H5Pset_fclose_degree(fapl, H5F_CLOSE_STRONG);
+ fh = H5Fopen(filename, H5F_ACC_RDONLY, fapl);
+ H5Pclose(fapl);
if ( fh < 0 ) {
ERROR("Couldn't open file (peaks/hdf5): %s\n", filename);
return NULL;
@@ -1486,6 +1538,7 @@ ImageFeatureList *image_hdf5_read_peaks_hdf5(const DataTemplate *dtempl,
if ( H5Sget_simple_extent_ndims(sh) != 2 ) {
ERROR("Peak list has the wrong dimensionality (%i).\n",
H5Sget_simple_extent_ndims(sh));
+ H5Sclose(sh);
close_hdf5(fh);
return NULL;
}
@@ -1860,7 +1913,7 @@ char **image_hdf5_expand_frames(const DataTemplate *dtempl,
{
char **path_evs;
int n_path_evs;
- hid_t fh;
+ hid_t fh, fapl;
int i;
int dims_expected;
struct ev_list full_evs;
@@ -1886,7 +1939,10 @@ char **image_hdf5_expand_frames(const DataTemplate *dtempl,
return NULL;
}
- fh = H5Fopen(filename, H5F_ACC_RDONLY, H5P_DEFAULT);
+ fapl = H5Pcreate(H5P_FILE_ACCESS);
+ H5Pset_fclose_degree(fapl, H5F_CLOSE_STRONG);
+ fh = H5Fopen(filename, H5F_ACC_RDONLY, fapl);
+ H5Pclose(fapl);
if ( fh < 0 ) {
ERROR("Couldn't open file (hdf5_expand_frames): %s\n", filename);
return NULL;
@@ -1950,6 +2006,7 @@ char **image_hdf5_expand_frames(const DataTemplate *dtempl,
"(%s has %i, expected %i)\n",
path, dims, dims_expected);
close_hdf5(fh);
+ H5Sclose(sh);
return NULL;
}
@@ -1958,15 +2015,19 @@ char **image_hdf5_expand_frames(const DataTemplate *dtempl,
if ( (size == NULL) || (placeholder_sizes == NULL) ) {
ERROR("Failed to allocate dimensions\n");
close_hdf5(fh);
+ H5Sclose(sh);
return NULL;
}
if ( H5Sget_simple_extent_dims(sh, size, NULL) < 0 ) {
ERROR("Failed to get size\n");
close_hdf5(fh);
+ H5Sclose(sh);
return NULL;
}
+ H5Sclose(sh);
+
n_placeholder_dims = 0;
for ( j=0; j<dims; j++ ) {
if ( p->dims[j] == DIM_PLACEHOLDER ) {
diff --git a/libcrystfel/src/image.c b/libcrystfel/src/image.c
index e91c09d6..3b215c08 100644
--- a/libcrystfel/src/image.c
+++ b/libcrystfel/src/image.c
@@ -1431,6 +1431,9 @@ void image_free(struct image *image)
for ( i=0; i<image->n_cached_headers; i++ ) {
cffree(image->header_cache[i]->header_name);
+ if ( image->header_cache[i]->type == HEADER_STR ) {
+ cffree(image->header_cache[i]->val_str);
+ }
cffree(image->header_cache[i]);
}
diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c
index 12cd9190..5b0f742f 100644
--- a/libcrystfel/src/index.c
+++ b/libcrystfel/src/index.c
@@ -50,6 +50,7 @@
#include "cell-utils.h"
#include "predict-refine.h"
#include "profile.h"
+#include "crystfel-mille.h"
#include "indexers/dirax.h"
#include "indexers/asdf.h"
#include "indexers/mosflm.h"
@@ -774,6 +775,12 @@ static int try_indexer(struct image *image, IndexingMethod indm,
}
profile_end("cell-compare-to-others");
+ if ( (mille != NULL) && !crystal_get_user_flag(cr) ) {
+ profile_start("mille-write");
+ crystfel_mille_write_record(mille);
+ profile_end("mille-write");
+ }
+
}
n_bad = remove_flagged_crystals(image);
diff --git a/libcrystfel/src/indexers/asdf.c b/libcrystfel/src/indexers/asdf.c
index c136efc1..6e1e5781 100644
--- a/libcrystfel/src/indexers/asdf.c
+++ b/libcrystfel/src/indexers/asdf.c
@@ -1308,7 +1308,7 @@ int run_asdf(struct image *image, void *ipriv)
}
-void *asdf_prepare(IndexingMethod *indm, UnitCell *cell)
+void *asdf_prepare(IndexingMethod *indm, UnitCell *cell, struct asdf_options *asdf_opts)
{
ERROR("This copy of CrystFEL was compiled without FFTW support.\n");
ERROR("To use asdf indexing, recompile with FFTW.\n");
diff --git a/libcrystfel/src/indexers/pinkindexer.c b/libcrystfel/src/indexers/pinkindexer.c
index 180246c2..14e18548 100644
--- a/libcrystfel/src/indexers/pinkindexer.c
+++ b/libcrystfel/src/indexers/pinkindexer.c
@@ -33,7 +33,6 @@
#include <stdlib.h>
-#include <sys/errno.h>
#include <argp.h>
#include "utils.h"
diff --git a/libcrystfel/src/indexers/taketwo.c b/libcrystfel/src/indexers/taketwo.c
index c87d6b3d..02fa9a23 100644
--- a/libcrystfel/src/indexers/taketwo.c
+++ b/libcrystfel/src/indexers/taketwo.c
@@ -522,84 +522,6 @@ static double matrix_trace(gsl_matrix *a)
return tr;
}
-static char *add_unique_axis(const char *inp, char ua)
-{
- char *pg = cfmalloc(64);
- if ( pg == NULL ) return NULL;
- snprintf(pg, 63, "%s_ua%c", inp, ua);
- return pg;
-}
-
-
-static char *get_chiral_holohedry(UnitCell *cell)
-{
- LatticeType lattice = cell_get_lattice_type(cell);
- char *pg;
- int add_ua = 1;
-
- switch (lattice)
- {
- case L_TRICLINIC:
- pg = "1";
- add_ua = 0;
- break;
-
- case L_MONOCLINIC:
- pg = "2";
- break;
-
- case L_ORTHORHOMBIC:
- pg = "222";
- add_ua = 0;
- break;
-
- case L_TETRAGONAL:
- pg = "422";
- break;
-
- case L_RHOMBOHEDRAL:
- pg = "3_R";
- add_ua = 0;
- break;
-
- case L_HEXAGONAL:
- if ( cell_get_centering(cell) == 'H' ) {
- pg = "3_H";
- add_ua = 0;
- } else {
- pg = "622";
- }
- break;
-
- case L_CUBIC:
- pg = "432";
- add_ua = 0;
- break;
-
- default:
- pg = "error";
- break;
- }
-
- if ( add_ua ) {
- return add_unique_axis(pg, cell_get_unique_axis(cell));
- } else {
- return cfstrdup(pg);
- }
-}
-
-
-static SymOpList *sym_ops_for_cell(UnitCell *cell)
-{
- SymOpList *rawList;
-
- char *pg = get_chiral_holohedry(cell);
- rawList = get_pointgroup(pg);
- cffree(pg);
-
- return rawList;
-}
-
static int rot_mats_are_similar(gsl_matrix *rot1, gsl_matrix *rot2,
gsl_matrix *sub, gsl_matrix *mul,
double *score, struct TakeTwoCell *cell)
@@ -1619,7 +1541,7 @@ static void set_gsl_matrix(gsl_matrix *mat,
static int generate_rotation_sym_ops(struct TakeTwoCell *ttCell)
{
- SymOpList *rawList = sym_ops_for_cell(ttCell->cell);
+ SymOpList *rawList = get_lattice_symmetry(ttCell->cell);
/* Now we must convert these into rotation matrices rather than hkl
* transformations (affects triclinic, monoclinic, rhombohedral and
@@ -1858,7 +1780,7 @@ static int gen_theoretical_vecs(UnitCell *cell, struct TheoryVec **cell_vecs,
&bsx, &bsy, &bsz,
&csx, &csy, &csz);
- SymOpList *rawList = sym_ops_for_cell(cell);
+ SymOpList *rawList = get_lattice_symmetry(cell);
cell_get_parameters(cell, &a, &b, &c, &alpha, &beta, &gamma);
diff --git a/libcrystfel/src/integration.c b/libcrystfel/src/integration.c
index d2c06db7..9ad22960 100644
--- a/libcrystfel/src/integration.c
+++ b/libcrystfel/src/integration.c
@@ -1467,11 +1467,6 @@ int integrate_rings_once(Reflection *refl,
* appropriate size */
if ( intensity < -aduph ) {
sig2_poisson = -aduph*intensity;
- } else if ( intensity < 0.0 ) {
- /* If the intensity is negative (by less than one
- * photon), assume the reflection is very weak and
- * therefore has a Poisson error of one photon. */
- sig2_poisson = aduph;
}
sigma = sqrt(sig2_poisson + bx->m*sig2_bg);
@@ -1816,6 +1811,8 @@ IntegrationMethod integration_method(const char *str, int *err)
meth = INTEGRATION_DEFAULTS_PROF2D;
} else if ( strcmp(methods[i], "none") == 0 ) {
+ cffree(methods[i]);
+ cffree(methods);
return INTEGRATION_NONE;
} else if ( strcmp(methods[i], "sat") == 0 ) {
diff --git a/libcrystfel/src/peakfinder8.c b/libcrystfel/src/peakfinder8.c
index 1def35d3..28bee8d7 100644
--- a/libcrystfel/src/peakfinder8.c
+++ b/libcrystfel/src/peakfinder8.c
@@ -141,7 +141,6 @@ static struct radial_stats_pixels *compute_rstats_pixels(struct radius_maps *rma
cffree(rsp);
return NULL;
}
- srand(0);
int n_pixels_per_bin = 100; // Can make this a parameter
diff --git a/libcrystfel/src/predict-refine.c b/libcrystfel/src/predict-refine.c
index 6da0d3d1..083d86c5 100644
--- a/libcrystfel/src/predict-refine.c
+++ b/libcrystfel/src/predict-refine.c
@@ -142,15 +142,15 @@ double r_gradient(int param, Reflection *refl, UnitCell *cell, double wavelength
/* Spot position gradients for diffraction physics (anything that changes the
* diffracted ray direction) */
int fs_ss_gradient_physics(int param, Reflection *refl, UnitCell *cell,
- struct detgeom_panel *p, gsl_matrix *Minv,
- double fs, double ss, double mu,
- float *fsg, float *ssg)
+ struct detgeom_panel *p, gsl_matrix *Minv,
+ double fs, double ss, double mu,
+ float *fsg, float *ssg)
{
signed int h, k, l;
gsl_vector *dRdp;
gsl_vector *v;
- get_indices(refl, &h, &k, &l);
+ get_symmetric_indices(refl, &h, &k, &l);
dRdp = gsl_vector_calloc(3);
@@ -302,7 +302,7 @@ int fs_ss_gradient(int param, Reflection *refl, UnitCell *cell,
double mu;
double fs, ss;
- get_indices(refl, &h, &k, &l);
+ get_symmetric_indices(refl, &h, &k, &l);
kpred = get_kpred(refl);
cell_get_reciprocal(cell, &asx, &asy, &asz,
&bsx, &bsy, &bsz,
@@ -811,7 +811,7 @@ int refine_prediction(struct image *image, Crystal *cr,
reflist = reflist_new();
n = pair_peaks(image, cr, reflist, rps);
- if ( n < 10 ) {
+ if ( n < 3 ) {
cffree(rps);
reflist_free(reflist);
return 1;
@@ -871,7 +871,8 @@ int refine_prediction(struct image *image, Crystal *cr,
total_shifts[0], total_shifts[1], total_shifts[2]);
crystal_add_notes(cr, tmp);
- if ( mille != NULL ) {
+ if ( (mille != NULL) && (n>4) ) {
+ crystfel_mille_delete_last_record(mille);
profile_start("mille-calc");
write_mille(mille, n, crystal_get_cell(cr), rps, image,
max_mille_level, Minvs);
@@ -887,18 +888,10 @@ int refine_prediction(struct image *image, Crystal *cr,
n = pair_peaks(image, cr, NULL, rps);
free_rps_noreflist(rps, n);
- if ( n < 10 ) {
- if ( mille != NULL ) {
- crystfel_mille_delete_last_record(mille);
- }
- return 1;
- }
- if ( mille != NULL ) {
- profile_start("mille-write");
- crystfel_mille_write_record(mille);
- profile_end("mille-write");
+ if ( n < 3 ) {
+ return 1;
+ } else {
+ return 0;
}
-
- return 0;
}
diff --git a/libcrystfel/src/utils.c b/libcrystfel/src/utils.c
index 7bfac487..a35ba272 100644
--- a/libcrystfel/src/utils.c
+++ b/libcrystfel/src/utils.c
@@ -135,39 +135,53 @@ gsl_matrix *matrix_invert(gsl_matrix *m)
perm = gsl_permutation_alloc(m->size1);
if ( perm == NULL ) {
ERROR("Couldn't allocate permutation\n");
- gsl_matrix_free(m);
return NULL;
}
inv = gsl_matrix_alloc(m->size1, m->size2);
if ( inv == NULL ) {
ERROR("Couldn't allocate inverse\n");
- gsl_matrix_free(m);
gsl_permutation_free(perm);
return NULL;
}
if ( gsl_linalg_LU_decomp(m, perm, &s) ) {
ERROR("Couldn't decompose matrix\n");
- gsl_matrix_free(m);
gsl_permutation_free(perm);
return NULL;
}
if ( gsl_linalg_LU_invert(m, perm, inv) ) {
- ERROR("Couldn't invert cell matrix:\n");
- gsl_matrix_free(m);
+ ERROR("Couldn't invert matrix:\n");
gsl_permutation_free(perm);
return NULL;
}
gsl_permutation_free(perm);
- gsl_matrix_free(m);
return inv;
}
+gsl_vector *solve_inv(gsl_vector *v, gsl_matrix *M)
+{
+ gsl_matrix *Minv;
+ gsl_vector *s;
+ int n;
+
+ n = v->size;
+ if ( v->size != M->size1 ) return NULL;
+ if ( v->size != M->size2 ) return NULL;
+
+ Minv = matrix_invert(M);
+ if ( Minv == NULL ) return NULL;
+
+ s = gsl_vector_calloc(n);
+ gsl_blas_dgemv(CblasNoTrans, 1.0, M, v, 0.0, s);
+ return s;
+}
+
+
static int check_eigen(gsl_vector *e_val, int verbose)
{
int i;
diff --git a/libcrystfel/src/utils.h b/libcrystfel/src/utils.h
index 30a7c4ef..08e059a7 100644
--- a/libcrystfel/src/utils.h
+++ b/libcrystfel/src/utils.h
@@ -75,6 +75,7 @@ extern "C" {
extern void show_matrix_eqn(gsl_matrix *M, gsl_vector *v);
extern void show_matrix(gsl_matrix *M);
extern void show_vector(gsl_vector *M);
+extern gsl_vector *solve_inv(gsl_vector *v, gsl_matrix *M);
extern gsl_vector *solve_svd(gsl_vector *v, gsl_matrix *M, int *n_filt,
int verbose);
extern gsl_matrix *matrix_mult2(gsl_matrix *A, gsl_matrix *B);