aboutsummaryrefslogtreecommitdiff
path: root/src/partialator.c
diff options
context:
space:
mode:
Diffstat (limited to 'src/partialator.c')
-rw-r--r--src/partialator.c177
1 files changed, 92 insertions, 85 deletions
diff --git a/src/partialator.c b/src/partialator.c
index 9b1803b0..7d3f19f6 100644
--- a/src/partialator.c
+++ b/src/partialator.c
@@ -169,14 +169,14 @@ static void add_to_csplit(struct custom_split *csplit, const char *id,
/* Write two-way split results (i.e. for CC1/2 etc) for this list of crystals */
-static void write_split(Crystal **crystals, int n_crystals, const char *outfile,
- int nthreads, PartialityModel pmodel,
+static void write_split(struct crystal_refls *crystals, int n_crystals,
+ const char *outfile, int nthreads, PartialityModel pmodel,
int min_measurements, SymOpList *sym, double push_res)
{
char tmp[1024];
RefList *split;
- Crystal **crystals1;
- Crystal **crystals2;
+ struct crystal_refls *crystals1;
+ struct crystal_refls *crystals2;
int n_crystals1 = 0;
int n_crystals2 = 0;
int i;
@@ -186,13 +186,12 @@ static void write_split(Crystal **crystals, int n_crystals, const char *outfile,
return;
}
- crystals1 = malloc(n_crystals * sizeof(Crystal *));
+ crystals1 = malloc(n_crystals * sizeof(struct crystal_refls));
if ( crystals1 == NULL ) return;
- crystals2 = malloc(n_crystals * sizeof(Crystal *));
+ crystals2 = malloc(n_crystals * sizeof(struct crystal_refls));
if ( crystals2 == NULL ) return;
-
for ( i=0; i<n_crystals; i++ ) {
if ( i % 2 ) {
crystals1[n_crystals1] = crystals[i];
@@ -259,14 +258,15 @@ static char *insert_into_filename(const char *fn, const char *add)
/* Write custom split results (including a two-way split) */
static void write_custom_split(struct custom_split *csplit, int dsn,
- Crystal **crystals, int n_crystals,
+ struct crystal_refls *crystals,
+ struct image **images, int n_crystals,
PartialityModel pmodel, int min_measurements,
double push_res, SymOpList *sym, int nthreads,
const char *outfile)
{
char *tmp;
RefList *split;
- Crystal *crystalsn[n_crystals];
+ struct crystal_refls crystalsn[n_crystals];
int n_crystalsn = 0;
int i;
@@ -277,8 +277,8 @@ static void write_custom_split(struct custom_split *csplit, int dsn,
char *id;
int dsn_crystal;
- fn = crystal_get_image(crystals[i])->filename;
- evs = crystal_get_image(crystals[i])->ev;
+ fn = images[i]->filename;
+ evs = images[i]->ev;
if ( evs == NULL ) evs = "//";
id = malloc(strlen(evs)+strlen(fn)+2);
@@ -362,7 +362,7 @@ static void show_help(const char *s)
}
-static signed int find_first_crystal(Crystal **crystals, int n_crystals,
+static signed int find_first_crystal(struct image **images, int n_crystals,
struct custom_split *csplit, int dsn)
{
int i;
@@ -373,8 +373,8 @@ static signed int find_first_crystal(Crystal **crystals, int n_crystals,
char *id;
int dsn_crystal;
- fn = crystal_get_image(crystals[i])->filename;
- evs = crystal_get_image(crystals[i])->ev;
+ fn = images[i]->filename;
+ evs = images[i]->ev;
if ( evs == NULL ) evs = "//";
id = malloc(strlen(evs)+strlen(fn)+2);
@@ -394,8 +394,8 @@ static signed int find_first_crystal(Crystal **crystals, int n_crystals,
}
-static void check_csplit(Crystal **crystals, int n_crystals,
- struct custom_split *csplit)
+static void check_csplit(struct crystal_refls *crystals, struct image **images,
+ int n_crystals, struct custom_split *csplit)
{
int i;
int n_nosplit = 0;
@@ -412,8 +412,8 @@ static void check_csplit(Crystal **crystals, int n_crystals,
char *id;
int dsn_crystal;
- fn = crystal_get_image(crystals[i])->filename;
- evs = crystal_get_image(crystals[i])->ev;
+ fn = images[i]->filename;
+ evs = images[i]->ev;
if ( evs == NULL ) evs = "//";
id = malloc(strlen(evs)+strlen(fn)+2);
@@ -437,7 +437,7 @@ static void check_csplit(Crystal **crystals, int n_crystals,
for ( i=0; i<csplit->n_datasets; i++ ) {
/* Try to find a crystal with dsn = i */
- if ( find_first_crystal(crystals, n_crystals, csplit, i) != -1 )
+ if ( find_first_crystal(images, n_crystals, csplit, i) != -1 )
{
n_cry++;
} else {
@@ -621,11 +621,10 @@ static void skip_to_end(FILE *fh)
}
-static int set_initial_params(Crystal *cr, FILE *fh, double force_bandwidth,
- double force_radius, double force_lambda)
+static int set_initial_params(Crystal *cr, struct image *image, FILE *fh,
+ double force_bandwidth, double force_radius,
+ double force_lambda)
{
- struct image *image = crystal_get_image(cr);
-
if ( fh != NULL ) {
int err;
@@ -736,7 +735,7 @@ static void write_to_pgraph(FILE *fh, RefList *list, RefList *full, Crystal *cr,
}
-static void write_pgraph(RefList *full, Crystal **crystals, int n_crystals,
+static void write_pgraph(RefList *full, struct crystal_refls *crystals, int n_crystals,
signed int iter, const char *suff,
const char *log_folder)
{
@@ -762,16 +761,16 @@ static void write_pgraph(RefList *full, Crystal **crystals, int n_crystals,
}
for ( i=0; i<n_crystals; i++ ) {
- if ( crystal_get_user_flag(crystals[i]) != 0 ) continue;
- write_to_pgraph(fh, crystal_get_reflections(crystals[i]), full,
- crystals[i], i, iter);
+ if ( crystal_get_user_flag(crystals[i].cr) != 0 ) continue;
+ write_to_pgraph(fh, crystals[i].refls, full,
+ crystals[i].cr, i, iter);
}
fclose(fh);
}
-static void all_residuals(Crystal **crystals, int n_crystals, RefList *full,
+static void all_residuals(struct crystal_refls *crystals, int n_crystals, RefList *full,
int no_free,
double *presidual, double *pfree_residual,
double *plog_residual, double *pfree_log_residual,
@@ -798,28 +797,28 @@ static void all_residuals(Crystal **crystals, int n_crystals, RefList *full,
double r, free_r, log_r, free_log_r;
int n;
- if ( crystal_get_user_flag(crystals[i]) ) continue;
+ if ( crystal_get_user_flag(crystals[i].cr) ) continue;
/* Scaling should have been done right before calling this */
- r = residual(crystals[i], full, 0, &n, NULL);
+ r = residual(crystals[i].refls, crystals[i].cr, full, 0, &n, NULL);
if ( n == 0 ) {
n_non_linear++;
} else if ( isnan(r) ) {
n_nan_linear++;
}
- free_r = residual(crystals[i], full, 1, &n, NULL);
+ free_r = residual(crystals[i].refls, crystals[i].cr, full, 1, &n, NULL);
if ( n == 0 ) {
n_non_linear_free++;
} else if ( isnan(free_r) ) {
n_nan_linear_free++;
}
- log_r = log_residual(crystals[i], full, 0, &n, NULL);
+ log_r = log_residual(crystals[i].refls, crystals[i].cr, full, 0, &n, NULL);
if ( n == 0 ) {
n_non_log++;
} else if ( isnan(log_r) ) {
n_nan_log++;
}
- free_log_r = log_residual(crystals[i], full, 1, &n, NULL);
+ free_log_r = log_residual(crystals[i].refls, crystals[i].cr, full, 1, &n, NULL);
if ( n == 0 ) {
n_non_log_free++;
} else if ( isnan(free_log_r) ) {
@@ -876,7 +875,7 @@ static void all_residuals(Crystal **crystals, int n_crystals, RefList *full,
}
-static void show_all_residuals(Crystal **crystals, int n_crystals,
+static void show_all_residuals(struct crystal_refls *crystals, int n_crystals,
RefList *full, int no_free)
{
double dev, free_dev, log_dev, free_log_dev;
@@ -896,7 +895,8 @@ struct log_qargs
{
int iter;
int next;
- Crystal **crystals;
+ struct crystal_refls *crystals;
+ struct image **images;
int n_crystals;
RefList *full;
int scaleflags;
@@ -909,6 +909,8 @@ struct log_qargs
struct log_args
{
Crystal *cr;
+ RefList *refls;
+ struct image *image;
RefList *full;
int scaleflags;
PartialityModel pmodel;
@@ -928,7 +930,9 @@ static void *get_log_task(void *vp)
task = malloc(sizeof(struct log_args));
if ( task == NULL ) return NULL;
- task->cr = qargs->crystals[qargs->next];
+ task->cr = qargs->crystals[qargs->next].cr;
+ task->refls = qargs->crystals[qargs->next].refls;
+ task->image = qargs->images[qargs->next];
task->full = qargs->full;
task->iter = qargs->iter;
task->cnum = qargs->next;
@@ -944,12 +948,13 @@ static void *get_log_task(void *vp)
static void write_logs(void *vp, int cookie)
{
struct log_args *args = vp;
- write_specgraph(args->cr, args->full, args->iter, args->cnum,
- args->log_folder);
- write_gridscan(args->cr, args->full, args->iter, args->cnum,
- args->scaleflags, args->pmodel, args->log_folder);
- write_test_logs(args->cr, args->full, args->iter, args->cnum,
- args->log_folder);
+ write_specgraph(args->refls, args->cr, args->image, args->full,
+ args->iter, args->cnum, args->log_folder);
+ write_gridscan(args->refls, args->cr, args->image, args->full,
+ args->iter, args->cnum, args->scaleflags, args->pmodel,
+ args->log_folder);
+ write_test_logs(args->cr, args->image, args->full, args->iter,
+ args->cnum, args->log_folder);
}
@@ -963,7 +968,8 @@ static void done_log(void *vqargs, void *vp)
}
-static void write_logs_parallel(Crystal **crystals, int n_crystals,
+static void write_logs_parallel(struct crystal_refls *crystals, struct image **images,
+ int n_crystals,
RefList *full, int iter, int n_threads,
int scaleflags, PartialityModel pmodel,
const char *log_folder)
@@ -974,6 +980,7 @@ static void write_logs_parallel(Crystal **crystals, int n_crystals,
qargs.next = 0;
qargs.full = full;
qargs.crystals = crystals;
+ qargs.images = images;
qargs.n_done = 0;
qargs.n_crystals = n_crystals;
qargs.scaleflags = scaleflags;
@@ -1084,7 +1091,8 @@ int main(int argc, char *argv[])
int no_scale = 0;
int no_Bscale = 0;
int no_pr = 0;
- Crystal **crystals;
+ struct crystal_refls *crystals;
+ struct image **images;
char *pmodel_str = NULL;
PartialityModel pmodel = PMODEL_XSPHERE;
int min_measurements = 2;
@@ -1526,6 +1534,7 @@ int main(int argc, char *argv[])
n_crystals = 0;
n_crystals_seen = 0;
crystals = NULL;
+ images = NULL;
if ( sparams_fn != NULL ) {
char line[1024];
sparams_fh = fopen(sparams_fn, "r");
@@ -1560,7 +1569,6 @@ int main(int argc, char *argv[])
do {
struct image *image;
- RefList *as;
int i;
image = stream_read_chunk(st, STREAM_REFLECTIONS);
@@ -1574,37 +1582,44 @@ int main(int argc, char *argv[])
for ( i=0; i<image->n_crystals; i++ ) {
Crystal *cr;
- Crystal **crystals_new;
+ struct crystal_refls *crystals_new;
+ struct image **images_new;
RefList *cr_refl;
- RefList *cr_refl_raw;
struct image *image_for_crystal;
double lowest_r;
n_crystals_seen++;
if ( n_crystals_seen <= start_after ) continue;
- if ( crystal_get_resolution_limit(image->crystals[i]) < min_res ) continue;
+ if ( crystal_get_resolution_limit(image->crystals[i].cr) < min_res ) continue;
- lowest_r = lowest_reflection(crystal_get_cell(image->crystals[i]));
- if ( crystal_get_profile_radius(image->crystals[i]) > 0.5*lowest_r ) {
+ lowest_r = lowest_reflection(crystal_get_cell(image->crystals[i].cr));
+ if ( crystal_get_profile_radius(image->crystals[i].cr) > 0.5*lowest_r ) {
ERROR("Rejecting %s %s crystal %i because "
"profile radius is obviously too big (%e %e).\n",
image->filename, image->ev, i,
- crystal_get_profile_radius(image->crystals[i]),
+ crystal_get_profile_radius(image->crystals[i].cr),
lowest_r);
continue;
}
crystals_new = realloc(crystals,
- (n_crystals+1)*sizeof(Crystal *));
+ (n_crystals+1)*sizeof(struct crystal_refls));
if ( crystals_new == NULL ) {
ERROR("Failed to allocate memory for crystal "
"list.\n");
return 1;
}
crystals = crystals_new;
- crystals[n_crystals] = crystal_copy_deep(image->crystals[i]);
- cr = crystals[n_crystals];
+ cr = crystal_copy(image->crystals[i].cr);
+ crystals[n_crystals].cr = cr;
+
+ images_new = realloc(images, (n_crystals+1)*sizeof(struct image *));
+ if ( images_new == NULL ) {
+ ERROR("Failed to allocate memory for image list\n");
+ return 1;
+ }
+ images = images_new;
/* Create a completely new, separate image
* structure for this crystal. */
@@ -1614,11 +1629,10 @@ int main(int argc, char *argv[])
return 1;
}
- crystal_set_image(cr, image_for_crystal);
*image_for_crystal = *image;
- image_for_crystal->n_crystals = 1;
- image_for_crystal->crystals = malloc(sizeof(Crystal *));
- image_for_crystal->crystals[0] = cr;
+ images[n_crystals] = image_for_crystal;
+ image_for_crystal->n_crystals = 0;
+ image_for_crystal->crystals = NULL;
image_for_crystal->filename = strdup(image->filename);
image_for_crystal->ev = safe_strdup(image->ev);
image_for_crystal->detgeom = NULL;
@@ -1629,20 +1643,18 @@ int main(int argc, char *argv[])
image_for_crystal->bad = NULL;
image_for_crystal->sat = NULL;
- /* This is the raw list of reflections */
- cr_refl_raw = crystal_get_reflections(cr);
-
- cr_refl = apply_max_adu(cr_refl_raw, max_adu);
- reflist_free(cr_refl_raw);
-
+ cr_refl = apply_max_adu(image->crystals[i].refls, max_adu);
if ( !no_free ) select_free_reflections(cr_refl, rng);
-
- as = asymmetric_indices(cr_refl, sym);
- crystal_set_reflections(cr, as);
- crystal_set_user_flag(cr, PRFLAG_OK);
+ image_add_crystal_refls(image_for_crystal, cr,
+ asymmetric_indices(cr_refl, sym));
reflist_free(cr_refl);
- if ( set_initial_params(cr, sparams_fh, force_bandwidth,
+ crystals[n_crystals].cr = image_for_crystal->crystals[0].cr;
+ crystals[n_crystals].refls = image_for_crystal->crystals[0].refls;
+
+ crystal_set_user_flag(cr, PRFLAG_OK);
+ if ( set_initial_params(cr, image_for_crystal,
+ sparams_fh, force_bandwidth,
force_radius, force_lambda) )
{
ERROR("Failed to set initial parameters\n");
@@ -1681,17 +1693,16 @@ int main(int argc, char *argv[])
STATUS("Initial partiality calculation...\n");
for ( icryst=0; icryst<n_crystals; icryst++ ) {
- Crystal *cr = crystals[icryst];
- update_predictions(cr);
+ Crystal *cr = crystals[icryst].cr;
+ RefList *refls = crystals[icryst].refls;
+ update_predictions(refls, cr, images[icryst]);
/* Polarisation correction requires kpred values */
- polarisation_correction(crystal_get_reflections(cr),
- crystal_get_cell(cr), polarisation);
-
- calculate_partialities(cr, pmodel);
+ polarisation_correction(refls, crystal_get_cell(cr), polarisation);
+ calculate_partialities(refls, cr, images[icryst], pmodel);
}
- if (csplit != NULL) check_csplit(crystals, n_crystals, csplit);
+ if (csplit != NULL) check_csplit(crystals, images, n_crystals, csplit);
/* Make a first pass at cutting out crap */
//STATUS("Early rejection...\n");
@@ -1716,7 +1727,7 @@ int main(int argc, char *argv[])
if ( do_write_logs ) {
write_pgraph(full, crystals, n_crystals, 0, "", log_folder);
- write_logs_parallel(crystals, n_crystals, full, 0, nthreads,
+ write_logs_parallel(crystals, images, n_crystals, full, 0, nthreads,
scaleflags, pmodel, log_folder);
}
@@ -1726,7 +1737,7 @@ int main(int argc, char *argv[])
STATUS("Scaling and refinement cycle %i of %i\n", itn+1, n_iter);
if ( !no_pr ) {
- refine_all(crystals, n_crystals, full, nthreads, pmodel,
+ refine_all(crystals, images, n_crystals, full, nthreads, pmodel,
itn+1, no_logs, sym, amb, scaleflags,
log_folder);
}
@@ -1771,7 +1782,7 @@ int main(int argc, char *argv[])
int j;
for ( j=0; j<csplit->n_datasets; j++ ) {
write_custom_split(csplit, j, crystals,
- n_crystals, pmodel,
+ images, n_crystals, pmodel,
min_measurements,
push_res, sym,
nthreads, tmp);
@@ -1805,7 +1816,7 @@ int main(int argc, char *argv[])
show_all_residuals(crystals, n_crystals, full, no_free);
if ( do_write_logs ) {
write_pgraph(full, crystals, n_crystals, -1, "", log_folder);
- write_logs_parallel(crystals, n_crystals, full, -1, nthreads,
+ write_logs_parallel(crystals, images, n_crystals, full, -1, nthreads,
scaleflags, pmodel, log_folder);
}
@@ -1827,7 +1838,7 @@ int main(int argc, char *argv[])
if ( csplit != NULL ) {
int i;
for ( i=0; i<csplit->n_datasets; i++ ) {
- write_custom_split(csplit, i, crystals, n_crystals,
+ write_custom_split(csplit, i, crystals, images, n_crystals,
pmodel, min_measurements, push_res,
sym, nthreads, outfile);
}
@@ -1835,10 +1846,6 @@ int main(int argc, char *argv[])
/* Clean up */
gsl_rng_free(rng);
- for ( icryst=0; icryst<n_crystals; icryst++ ) {
- struct image *image = crystal_get_image(crystals[icryst]);
- image_free(image);
- }
free_contribs(full);
reflist_free(full);
free_symoplist(sym);