Commit 4ca2f404 authored by Thomas White's avatar Thomas White Committed by Thomas White
Browse files

Implement "multi-indexing"

parent cc323ec8
......@@ -43,39 +43,67 @@ static IndexingPrivate *indexing_private(IndexingMethod indm)
}
IndexingPrivate *prepare_indexing(IndexingMethod indm, UnitCell *cell,
IndexingPrivate **prepare_indexing(IndexingMethod *indm, UnitCell *cell,
const char *filename, struct detector *det,
double nominal_photon_energy)
{
switch ( indm ) {
case INDEXING_NONE :
return indexing_private(indm);
case INDEXING_DIRAX :
return indexing_private(indm);
case INDEXING_MOSFLM :
return indexing_private(indm);
case INDEXING_TEMPLATE :
return generate_templates(cell, filename, det,
nominal_photon_energy);
int n;
int nm = 0;
IndexingPrivate **iprivs;
while ( indm[nm] != INDEXING_NONE ) nm++;
STATUS("Preparing %i indexing methods.\n", nm);
iprivs = malloc((nm+1) * sizeof(IndexingPrivate *));
for ( n=0; n<nm; n++ ) {
switch ( indm[n] ) {
case INDEXING_NONE :
ERROR("Tried to prepare INDEXING_NONE!\n");
break;
case INDEXING_DIRAX :
iprivs[n] = indexing_private(indm[n]);
break;
case INDEXING_MOSFLM :
iprivs[n] = indexing_private(indm[n]);
break;
case INDEXING_TEMPLATE :
iprivs[n] = generate_templates(cell, filename, det,
nominal_photon_energy);
break;
}
n++;
}
return 0;
iprivs[n] = NULL;
return iprivs;
}
void cleanup_indexing(IndexingPrivate *priv)
void cleanup_indexing(IndexingPrivate **priv)
{
switch ( priv->indm ) {
case INDEXING_NONE :
free(priv);
break;
case INDEXING_DIRAX :
free(priv);
break;
case INDEXING_MOSFLM :
free(priv);
break;
case INDEXING_TEMPLATE :
free_templates(priv);
int n = 0;
while ( priv[n] != NULL ) {
switch ( priv[n]->indm ) {
case INDEXING_NONE :
free(priv[n]);
break;
case INDEXING_DIRAX :
free(priv[n]);
break;
case INDEXING_MOSFLM :
free(priv[n]);
break;
case INDEXING_TEMPLATE :
free_templates(priv[n]);
}
n++;
}
}
......@@ -246,81 +274,93 @@ void map_all_peaks(struct image *image)
}
void index_pattern(struct image *image, UnitCell *cell, IndexingMethod indm,
int cellr, int verbose, IndexingPrivate *ipriv)
void index_pattern(struct image *image, UnitCell *cell, IndexingMethod *indm,
int cellr, int verbose, IndexingPrivate **ipriv)
{
int i;
int n = 0;
map_all_peaks(image);
write_drx(image);
image->ncells = 0;
while ( indm[n] != INDEXING_NONE ) {
/* Index (or not) as appropriate */
switch ( indm ) {
case INDEXING_NONE :
return;
case INDEXING_DIRAX :
run_dirax(image);
break;
case INDEXING_MOSFLM :
write_spt(image);
write_img(image); /* dummy image. not needed for each frame.*/
run_mosflm(image,cell);
break;
case INDEXING_TEMPLATE :
match_templates(image, ipriv);
break;
}
if ( image->ncells == 0 ) {
STATUS("No candidate cells found.\n");
return;
}
image->ncells = 0;
if ( (cellr == CELLR_NONE) || (indm == INDEXING_TEMPLATE) ) {
image->indexed_cell = image->candidate_cells[0];
if ( verbose ) {
STATUS("--------------------\n");
STATUS("The indexed cell (matching not performed):\n");
cell_print(image->indexed_cell);
STATUS("--------------------\n");
/* Index as appropriate */
switch ( indm[n] ) {
case INDEXING_NONE :
return;
case INDEXING_DIRAX :
STATUS("Running DirAx...\n");
write_drx(image);
run_dirax(image);
break;
case INDEXING_MOSFLM :
STATUS("Running MOSFLM...\n");
write_spt(image);
write_img(image); /* dummy image */
run_mosflm(image, cell);
break;
case INDEXING_TEMPLATE :
match_templates(image, ipriv[n]);
break;
}
return;
}
for ( i=0; i<image->ncells; i++ ) {
UnitCell *new_cell = NULL;
if ( verbose ) {
STATUS("--------------------\n");
STATUS("Candidate cell %i (before matching):\n", i);
cell_print(image->candidate_cells[i]);
STATUS("--------------------\n");
if ( image->ncells == 0 ) {
STATUS("No candidate cells found.\n");
return;
}
/* Match or reduce the cell as appropriate */
switch ( cellr ) {
case CELLR_NONE :
/* Never happens */
break;
case CELLR_REDUCE :
new_cell = match_cell(image->candidate_cells[i],
cell, verbose, 1);
break;
case CELLR_COMPARE :
new_cell = match_cell(image->candidate_cells[i],
cell, verbose, 0);
break;
if ( (cellr == CELLR_NONE) || (indm[n] == INDEXING_TEMPLATE) ) {
image->indexed_cell = image->candidate_cells[0];
if ( verbose ) {
STATUS("--------------------\n");
STATUS("The indexed cell (matching not"
" performed):\n");
cell_print(image->indexed_cell);
STATUS("--------------------\n");
}
return;
}
image->indexed_cell = new_cell;
if ( new_cell != NULL ) {
STATUS("Matched on attempt %i.\n", i);
goto done;
for ( i=0; i<image->ncells; i++ ) {
UnitCell *new_cell = NULL;
if ( verbose ) {
STATUS("--------------------\n");
STATUS("Candidate cell %i (before matching):\n",
i);
cell_print(image->candidate_cells[i]);
STATUS("--------------------\n");
}
/* Match or reduce the cell as appropriate */
switch ( cellr ) {
case CELLR_NONE :
/* Never happens */
break;
case CELLR_REDUCE :
new_cell = match_cell(image->candidate_cells[i],
cell, verbose, 1);
break;
case CELLR_COMPARE :
new_cell = match_cell(image->candidate_cells[i],
cell, verbose, 0);
break;
}
image->indexed_cell = new_cell;
if ( new_cell != NULL ) {
STATUS("Matched on attempt %i.\n", i);
goto done;
}
}
/* Move on to the next indexing method */
n++;
}
done:
......@@ -328,3 +368,38 @@ done:
cell_free(image->candidate_cells[i]);
}
}
IndexingMethod *build_indexer_list(const char *str, int *need_cell)
{
int n, i;
char **methods;
IndexingMethod *list;
*need_cell = 0;
n = assplode(str, ",", &methods, ASSPLODE_NONE);
list = malloc((n+1)*sizeof(IndexingMethod));
for ( i=0; i<n; i++ ) {
if ( strcmp(methods[i], "dirax") == 0) {
list[i] = INDEXING_DIRAX;
} else if ( strcmp(methods[i], "mosflm") == 0) {
list[i] = INDEXING_MOSFLM;
} else if ( strcmp(methods[i], "template") == 0) {
list[i] = INDEXING_TEMPLATE;
*need_cell = 1;
} else {
ERROR("Unrecognised indexing method '%s'\n",
methods[i]);
return NULL;
}
free(methods[i]);
}
free(methods);
list[i] = INDEXING_NONE;
return list;
}
......@@ -42,7 +42,7 @@ enum {
typedef struct _indexingprivate IndexingPrivate;
extern IndexingPrivate *prepare_indexing(IndexingMethod indm, UnitCell *cell,
extern IndexingPrivate **prepare_indexing(IndexingMethod *indm, UnitCell *cell,
const char *filename,
struct detector *det,
double nominal_photon_energy);
......@@ -50,9 +50,11 @@ extern IndexingPrivate *prepare_indexing(IndexingMethod indm, UnitCell *cell,
extern void map_all_peaks(struct image *image);
extern void index_pattern(struct image *image, UnitCell *cell,
IndexingMethod indm, int cellr, int verbose,
IndexingPrivate *priv);
IndexingMethod *indm, int cellr, int verbose,
IndexingPrivate **priv);
extern void cleanup_indexing(IndexingPrivate *priv);
extern void cleanup_indexing(IndexingPrivate **priv);
extern IndexingMethod *build_indexer_list(const char *str, int *need_cell);
#endif /* INDEX_H */
......@@ -68,8 +68,8 @@ struct static_index_args
float threshold;
float min_gradient;
struct detector *det;
IndexingMethod indm;
IndexingPrivate *ipriv;
IndexingMethod *indm;
IndexingPrivate **ipriv;
const double *intensities;
const unsigned char *flags;
const char *sym; /* Symmetry of "intensities" and "flags" */
......@@ -325,7 +325,7 @@ static void process_image(void *pp, int cookie)
int config_gpu = pargs->static_args.config_gpu;
int config_simulate = pargs->static_args.config_simulate;
int config_polar = pargs->static_args.config_polar;
IndexingMethod indm = pargs->static_args.indm;
IndexingMethod *indm = pargs->static_args.indm;
const double *intensities = pargs->static_args.intensities;
const unsigned char *flags = pargs->static_args.flags;
struct gpu_context *gctx = pargs->static_args.gctx;
......@@ -543,7 +543,10 @@ int main(int argc, char *argv[])
float min_gradient = 100000.0;
struct detector *det;
char *geometry = NULL;
IndexingMethod indm;
IndexingMethod *indm;
IndexingPrivate **ipriv;
int indexer_needs_cell;
int reduction_needs_cell;
char *indm_str = NULL;
UnitCell *cell;
double *intensities = NULL;
......@@ -561,7 +564,6 @@ int main(int argc, char *argv[])
pthread_mutex_t gpu_mutex = PTHREAD_MUTEX_INITIALIZER;
char prepare_line[1024];
char prepare_filename[1024];
IndexingPrivate *ipriv;
struct queue_args qargs;
struct beam_params *beam = NULL;
double nominal_photon_energy;
......@@ -784,31 +786,30 @@ int main(int argc, char *argv[])
" try to index anything.\n"
"If that isn't what you wanted, re-run with"
" --indexing=<method>.\n");
indm = INDEXING_NONE;
} else if ( strcmp(indm_str, "none") == 0 ) {
indm = INDEXING_NONE;
} else if ( strcmp(indm_str, "dirax") == 0) {
indm = INDEXING_DIRAX;
} else if ( strcmp(indm_str, "mosflm") == 0) {
indm = INDEXING_MOSFLM;
} else if ( strcmp(indm_str, "template") == 0) {
indm = INDEXING_TEMPLATE;
indm = NULL;
} else {
ERROR("Unrecognised indexing method '%s'\n", indm_str);
return 1;
indm = build_indexer_list(indm_str, &indexer_needs_cell);
if ( indm == NULL ) {
ERROR("Invalid indexer list '%s'\n", indm_str);
return 1;
}
free(indm_str);
}
free(indm_str);
reduction_needs_cell = 0;
if ( scellr == NULL ) {
STATUS("You didn't specify a cell reduction method, so I'm"
" going to use 'reduce'.\n");
cellr = CELLR_REDUCE;
reduction_needs_cell = 1;
} else if ( strcmp(scellr, "none") == 0 ) {
cellr = CELLR_NONE;
} else if ( strcmp(scellr, "reduce") == 0) {
cellr = CELLR_REDUCE;
reduction_needs_cell = 1;
} else if ( strcmp(scellr, "compare") == 0) {
cellr = CELLR_COMPARE;
reduction_needs_cell = 1;
} else {
ERROR("Unrecognised cell reduction method '%s'\n", scellr);
return 1;
......@@ -827,14 +828,15 @@ int main(int argc, char *argv[])
}
free(geometry);
if ( (cellr != CELLR_NONE) || (indm == INDEXING_TEMPLATE) ) {
if ( reduction_needs_cell || indexer_needs_cell ) {
cell = load_cell_from_pdb(pdb);
if ( cell == NULL ) {
ERROR("Couldn't read unit cell (from %s)\n", pdb);
return 1;
}
} else {
STATUS("No cell needed because --no-match was used.\n");
STATUS("No cell needed for these choices of indexing"
" and reduction.\n");
cell = NULL;
}
free(pdb);
......@@ -866,11 +868,15 @@ int main(int argc, char *argv[])
qargs.use_this_one_instead = prepare_line;
/* Prepare the indexer */
ipriv = prepare_indexing(indm, cell, prepare_filename, det,
nominal_photon_energy);
if ( ipriv == NULL ) {
ERROR("Failed to prepare indexing.\n");
return 1;
if ( indm != NULL ) {
ipriv = prepare_indexing(indm, cell, prepare_filename, det,
nominal_photon_energy);
if ( ipriv == NULL ) {
ERROR("Failed to prepare indexing.\n");
return 1;
}
} else {
ipriv = NULL;
}
gsl_set_error_handler_off();
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment