From 3926cfee958f2b3687997f22881c6d2593bbd3b7 Mon Sep 17 00:00:00 2001 From: Pascal Hogan-Lamarre Date: Mon, 14 Sep 2020 23:58:00 -0400 Subject: Removed global variable, include profile_radius and diff_limit in sol file, cleaned up warning messages, not working yet for multi-crystals --- libcrystfel/src/index.c | 265 ++++++++++++++++++++++++++---------------------- 1 file changed, 145 insertions(+), 120 deletions(-) diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c index 2b5907ba..0150787a 100644 --- a/libcrystfel/src/index.c +++ b/libcrystfel/src/index.c @@ -42,6 +42,7 @@ #include #include #include +#include #include "image.h" #include "utils.h" @@ -65,28 +66,28 @@ /** \file index.h */ +#define NPARAMS_PER_LINE 13 // There are 9 vector components, 2 detector shifts, 1 profile radius, 1 resolution limit +#define NKEYS_PER_LINE 4 // The keys are the filename, event path, event dim and crystal number + struct skip_private { - char path_to_sol[50]; - UnitCell *cellTemplate; - float solutions[]; + UnitCell *cellTemplate; + struct record_t *sol_hash; }; -struct record_key_t{ //added +struct record_key_t{ //Hash table keys char filename[100]; - char event_path[35]; + char event_path[100]; int event_dim; + int crystal_number; }; -struct record_t{ //added +struct record_t{ //Hash table struct record_key_t key; - int line; + float solution[NPARAMS_PER_LINE]; UT_hash_handle hh; }; -//declared as global variable -struct record_t *line_hash = NULL; //added - struct _indexingprivate { IndexingFlags flags; @@ -98,22 +99,24 @@ struct _indexingprivate void **engine_private; }; -void delete_all() { - struct record_t *current_line, *tmp; +void print_struct(struct record_t *sol_hash) { + struct record_t *s; + s = (struct record_t *)malloc(sizeof *s); + memset(s, 0, sizeof *s); - HASH_ITER(hh, line_hash, current_line, tmp) { - HASH_DEL(line_hash, current_line); - free(current_line); - } + for(s=sol_hash; s != NULL; s=(struct record_t*)(s->hh.next)) { + printf("This entry corresponds to event %s//%d and crystal_number %d in file %s \n", s->key.event_path, s->key.event_dim, s->key.crystal_number, s->key.filename); + } } -void print_struct() { +void full_print_struct(struct record_t *sol_hash) { struct record_t *s; s = (struct record_t *)malloc(sizeof *s); memset(s, 0, sizeof *s); - for(s=line_hash; s != NULL; s=(struct record_t*)(s->hh.next)) { - printf("The line is event_path %s and event_dim %d in file %s \n", s->key.event_path, s->key.event_dim, s->key.filename); + for(s=sol_hash; s != NULL; s=(struct record_t*)(s->hh.next)) { + printf("This entry corresponds to event %s//%d and crystal_number %d in file %s \n", s->key.event_path, s->key.event_dim, s->key.crystal_number, s->key.filename); + printf("The solution parameters are %e %e %e %e %e %e %e %e %e %e %e %e %e \n", s->solution[0], s->solution[1], s->solution[2], s->solution[3], s->solution[4], s->solution[5], s->solution[6], s->solution[7], s->solution[8], s->solution[9], s->solution[10], s->solution[11], s->solution[12]); } } @@ -137,7 +140,6 @@ static void set_last_task(char *lt, const char *task) strcpy(lt, task); } - static void show_indexing_flags(IndexingFlags flags) { STATUS("Indexing parameters:\n"); @@ -153,30 +155,60 @@ static void show_indexing_flags(IndexingFlags flags) onoff(flags & INDEXING_RETRY)); } +int ncrystals_in_sol(char *path) +{ + FILE *fh; + int count = 0; // Line counter (result) + char c; // To store a character read from file + + fh = fopen(path, "r"); + + if ( fh == NULL ) { + ERROR("%s not found by ncrystals_in_sol\n",path); + return 0; + } + + for (c = getc(fh); c != EOF; c = getc(fh)) + if (c == '\n') // Increment count if this character is newline + count = count + 1; + + //For the last line, which has no \n at the end + count = count + 1; + + // Close the file + fclose(fh); + + // STATUS("Found indexing file %s containing %d lines. \n", path, count); + + return count; + +} + void *skip_prepare(char *solution_filename, UnitCell *cell) { FILE *fh; int nlines; int nparams_in_solution; - int nparams_per_line; int nentries; - char *filename[35]; //??? is that correct - char *event_path[35]; //??? is that correct - int event_dim; + char filename[100]; + char event_path[100]; + int event_dim; + int crystal_number; int current_line; int position_in_current_line; - struct record_t *item; + struct record_t *sol_hash = NULL; + struct record_t *item = NULL; + float params[NPARAMS_PER_LINE]; + char path_to_sol[50],extension[10]; + char cwd[PATH_MAX]; //Assembling solution file name from input file name to crystfel - char path_to_sol[50],extension[10]; strcpy(path_to_sol,"../"); strcpy(extension,".sol"); - strtok(solution_filename, "."); //Takes the first part of the list name (before the .) strcat(path_to_sol, strtok(solution_filename, ".") ); //Add .sol at the previously splitted name strcat(path_to_sol, extension ); - char cwd[PATH_MAX]; if (getcwd(cwd, sizeof(cwd)) != NULL) { //printf("Path where skip_prepare is ran: %s\n", cwd); } @@ -192,22 +224,19 @@ void *skip_prepare(char *solution_filename, UnitCell *cell) } nlines = ncrystals_in_sol(path_to_sol); - nparams_per_line = 11; /* There are 9 vector components and 2 detector shifts */ - nparams_in_solution = nlines*nparams_per_line; /* total crystal parameters in solution file */ - nentries = nlines*(nparams_per_line+3); /* total entries in solution file */ + nparams_in_solution = nlines*NPARAMS_PER_LINE; /* total crystal parameters in solution file */ + nentries = nlines*(NPARAMS_PER_LINE+NKEYS_PER_LINE); /* total entries in solution file */ - float *params = malloc(nparams_in_solution * sizeof( float)); - STATUS("Parsing solution file containing %d lines...\n", nlines); //Reads indexing solutions - int j = 0; //index that follows the solution parameter + int j = 0; //index that follows the solution parameter [0,NPARAMS_PER_LINE] for(int i = 0; i < nentries; i++) { - current_line = i/(nparams_per_line+3); + current_line = i/(NPARAMS_PER_LINE+NKEYS_PER_LINE); - position_in_current_line = (i)%(nparams_per_line+3); + position_in_current_line = (i)%(NPARAMS_PER_LINE+NKEYS_PER_LINE); if (position_in_current_line == 0){ @@ -225,106 +254,84 @@ void *skip_prepare(char *solution_filename, UnitCell *cell) if (position_in_current_line == 1){ if ( fscanf(fh, "%s", event_path) != 1 ) { - printf("Failed to read an event_path from indexing.sol with skip_prepare\n"); + printf("Failed to read an event path from indexing.sol with skip_prepare\n"); return 0; } } - + if (position_in_current_line == 2){ + if ( fscanf(fh, "%d", &event_dim) != 1 ) { - printf("Failed to read an event from indexing.sol with skip_prepare\n"); + printf("Failed to read an event dim from indexing.sol with skip_prepare\n"); + return 0; + } + } + + if (position_in_current_line == 3){ + if ( fscanf(fh, "%d", &crystal_number) != 1 ) { + printf("Failed to read a crystal number from indexing.sol with skip_prepare\n"); + return 0; + } + } + + if (position_in_current_line > 3){ + if ( fscanf(fh, "%e", ¶ms[j]) != 1 ) { + printf("Failed to read a parameter from indexing.sol with skip_prepare\n"); return 0; } + j+=1; + + + } + + if (j == (NPARAMS_PER_LINE) ){ - //Prepare to add to the hash table once an event is read + //Prepare to add to the hash table item = (struct record_t *)malloc(sizeof *item); memset(item, 0, sizeof *item); strcpy(item->key.filename, filename); strcpy(item->key.event_path, event_path); item->key.event_dim = event_dim; - item->line=current_line; + item->key.crystal_number = crystal_number; + for (int k = 0; k < NPARAMS_PER_LINE; k++){ + item->solution[k] = params[k]; + } //Verify the uniqueness of the key struct record_t *uniqueness_test; - HASH_FIND(hh, line_hash, &item->key, sizeof(struct record_key_t), uniqueness_test); /* id already in the hash? */ + HASH_FIND(hh, sol_hash, &item->key, sizeof(struct record_key_t), uniqueness_test); // id already in the hash? if (uniqueness_test==NULL) { //If not in the table, add - HASH_ADD(hh, line_hash, key, sizeof(struct record_key_t), item); + HASH_ADD(hh, sol_hash, key, sizeof(struct record_key_t), item); } else{ //If already in the table, break with error - printf("Keys to the data must be unique! Verify the combination filename, event_path, event_dim"); + printf("Keys to the data must be unique! Verify the combination filename, event, crystal_number"); return 0; } + j=0; - - } - - if (position_in_current_line > 2){ - if ( fscanf(fh, "%e", ¶ms[j]) != 1 ) { - printf("Failed to read a parameter from indexing.sol with skip_prepare\n"); - return 0; - } - j+=1; } } fclose(fh); - + STATUS("Solution file parsing done. Have %d parameters and %d total entries.\n", nparams_in_solution, nentries) struct skip_private *dp; - dp = (struct skip_private *) malloc( sizeof(struct skip_private) + nentries * sizeof( float) ); + dp = (struct skip_private *) malloc( sizeof(struct skip_private)); if ( dp == NULL ) return NULL; - - memcpy(dp->path_to_sol, path_to_sol, sizeof path_to_sol); - - for (int k = 0; k < nparams_in_solution; k++) - { - dp->solutions[k] = params[k]; - } - dp->cellTemplate = cell; - - free(params); + dp->sol_hash = sol_hash; STATUS("Solution lookup table initialized!\n"); return (void *)dp; } -int ncrystals_in_sol(char *path) -{ - FILE *fh; - int count = 0; // Line counter (result) - char c; // To store a character read from file - - fh = fopen(path, "r"); - - if ( fh == NULL ) { - ERROR("%s not found by ncrystals_in_sol\n",path); - return 0; - } - - for (c = getc(fh); c != EOF; c = getc(fh)) - if (c == '\n') // Increment count if this character is newline - count = count + 1; - - //For the last line, which has no \n at the end - count = count + 1; - - // Close the file - fclose(fh); - - // STATUS("Found indexing file %s containing %d lines. \n", path, count); - - return count; - -} - //Fonction from pinkindexer.c to update the detector center of individual crystals //hack for electron crystallography while crystal_set_det_shift is not working approprietly static void update_detector(struct detector *det, double xoffs, double yoffs) @@ -338,15 +345,13 @@ static void update_detector(struct detector *det, double xoffs, double yoffs) } } -static int skip_index(struct image *image, void *mpriv) +static int skip_index(struct image *image, void *mpriv, int crystal_number) { Crystal *cr; UnitCell *cell; - float asx, asy, asz, bsx, bsy, bsz, csx, csy, csz, xshift, yshift; - int serial; - int nparams_per_line; - struct record_t *item, *p; - int line_number; + float asx, asy, asz, bsx, bsy, bsz, csx, csy, csz, xshift, yshift, profile_radius, resolution_limit; + struct record_t *item, *p, *pprime; + float *sol; struct skip_private *dp = (struct skip_private *)mpriv; @@ -356,29 +361,29 @@ static int skip_index(struct image *image, void *mpriv) strcpy(item->key.filename, image->filename); strcpy(item->key.event_path, *image->event->path_entries); item->key.event_dim = *image->event->dim_entries; + item->key.crystal_number = crystal_number; - HASH_FIND(hh, line_hash, &item->key, sizeof(struct record_key_t), p); /* id already in the hash? */ + HASH_FIND(hh, dp->sol_hash, &item->key, sizeof(struct record_key_t), p); //id already in the hash? if (p==NULL) { - /* - STATUS("Not indexing file %s, event %s %d \n", image->filename, *image->event->path_entries, *image->event->dim_entries); - */ + //STATUS("Not indexing file %s, event %s %d \n", image->filename, *image->event->path_entries, *image->event->dim_entries); return 0; } - line_number = p->line; - nparams_per_line = 11; - - asx = dp->solutions[(line_number * nparams_per_line) + 0]; - asy = dp->solutions[(line_number * nparams_per_line) + 1]; - asz = dp->solutions[(line_number * nparams_per_line) + 2]; - bsx = dp->solutions[(line_number * nparams_per_line) + 3]; - bsy = dp->solutions[(line_number * nparams_per_line) + 4]; - bsz = dp->solutions[(line_number * nparams_per_line) + 5]; - csx = dp->solutions[(line_number * nparams_per_line) + 6]; - csy = dp->solutions[(line_number * nparams_per_line) + 7]; - csz = dp->solutions[(line_number * nparams_per_line) + 8]; - xshift = dp->solutions[(line_number * nparams_per_line) + 9]; - yshift = dp->solutions[(line_number * nparams_per_line) + 10]; + sol = &(p->solution)[0]; + + asx = sol[0]; + asy = sol[1]; + asz = sol[2]; + bsx = sol[3]; + bsy = sol[4]; + bsz = sol[5]; + csx = sol[6]; + csy = sol[7]; + csz = sol[8]; + xshift = sol[9]; + yshift = sol[10]; + profile_radius = sol[11]; + resolution_limit = sol[12]; cr = crystal_new(); cell = cell_new(); @@ -390,9 +395,29 @@ static int skip_index(struct image *image, void *mpriv) crystal_set_cell(cr, cell); crystal_set_det_shift(cr, xshift * 1e-3, yshift * 1e-3); update_detector(image->det, xshift * 1e-3, yshift * 1e-3); + crystal_set_profile_radius(cr, profile_radius*1e9); + crystal_set_resolution_limit(cr, resolution_limit*1e9); image_add_crystal(image, cr); + //Look for additional crystals + item->key.crystal_number = crystal_number+1; + HASH_FIND(hh, dp->sol_hash, &item->key, sizeof(struct record_key_t), pprime); //id already in the hash? + if (pprime==NULL) { + + //If no more crystal, done + return 1; + } + else{ + //If more crystals, recursive call for next crystal in line + skip_index(image,mpriv,crystal_number+1); + } + + + dp->sol_hash = NULL; //Clean up local copy of the hash table + + return 1; + } @@ -771,7 +796,6 @@ void cleanup_indexing(IndexingPrivate *ipriv) case INDEXING_FILE : free(ipriv->engine_private[n]); - delete_all(); break; case INDEXING_TAKETWO : @@ -899,7 +923,8 @@ static int try_indexer(struct image *image, IndexingMethod indm, case INDEXING_FILE : set_last_task(last_task, "indexing:file"); - r = skip_index(image,mpriv); + int crystal_number = 0; + r = skip_index(image,mpriv,crystal_number); break; case INDEXING_FELIX : -- cgit v1.2.3