aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--libcrystfel/src/cell.c2
-rw-r--r--libcrystfel/src/cell.h2
-rw-r--r--libcrystfel/src/index.c2
-rw-r--r--libcrystfel/src/index.h3
-rw-r--r--src/indexamajig.c385
5 files changed, 225 insertions, 169 deletions
diff --git a/libcrystfel/src/cell.c b/libcrystfel/src/cell.c
index 8412fc62..2e87b770 100644
--- a/libcrystfel/src/cell.c
+++ b/libcrystfel/src/cell.c
@@ -690,7 +690,7 @@ static int same_vector(struct cvec a, struct cvec b)
/* Attempt to make 'cell' fit into 'template' somehow */
UnitCell *match_cell(UnitCell *cell, UnitCell *template, int verbose,
- float *tols, int reduce)
+ const float *tols, int reduce)
{
signed int n1l, n2l, n3l;
double asx, asy, asz;
diff --git a/libcrystfel/src/cell.h b/libcrystfel/src/cell.h
index cbe84772..bd2719dd 100644
--- a/libcrystfel/src/cell.h
+++ b/libcrystfel/src/cell.h
@@ -118,7 +118,7 @@ extern UnitCell *rotate_cell(UnitCell *in, double omega, double phi,
extern void cell_print(UnitCell *cell);
extern UnitCell *match_cell(UnitCell *cell, UnitCell *tempcell, int verbose,
- float *ltl, int reduce);
+ const float *ltl, int reduce);
extern UnitCell *match_cell_ab(UnitCell *cell, UnitCell *tempcell);
diff --git a/libcrystfel/src/index.c b/libcrystfel/src/index.c
index 3d0f164b..2baa2644 100644
--- a/libcrystfel/src/index.c
+++ b/libcrystfel/src/index.c
@@ -157,7 +157,7 @@ void map_all_peaks(struct image *image)
void index_pattern(struct image *image, UnitCell *cell, IndexingMethod *indm,
int cellr, int verbose, IndexingPrivate **ipriv,
- int config_insane, float *ltl)
+ int config_insane, const float *ltl)
{
int i;
int n = 0;
diff --git a/libcrystfel/src/index.h b/libcrystfel/src/index.h
index c4e15f09..9d23f3fb 100644
--- a/libcrystfel/src/index.h
+++ b/libcrystfel/src/index.h
@@ -74,7 +74,8 @@ 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, int config_insane, float *ltl);
+ IndexingPrivate **priv, int config_insane,
+ const float *ltl);
extern void cleanup_indexing(IndexingPrivate **priv);
diff --git a/src/indexamajig.c b/src/indexamajig.c
index ed84c6a8..75fcb1f5 100644
--- a/src/indexamajig.c
+++ b/src/indexamajig.c
@@ -70,10 +70,6 @@
/* Write statistics at APPROXIMATELY this interval */
#define STATS_EVERY_N_SECONDS (5)
-#define LINE_LENGTH 1024
-
-#define BUFFER PIPE_BUF
-
enum {
PEAK_ZAEF,
PEAK_HDF5,
@@ -228,72 +224,73 @@ static void show_help(const char *s)
static char *get_pattern(FILE *fh)
{
char *rval;
- char line[LINE_LENGTH];
- rval = fgets(line, LINE_LENGTH - 1, fh);
+ char *line;
+
+ line = malloc(1024);
+ if ( line == NULL ) {
+ ERROR("Couldn't allocate memory for filename\n");
+ return NULL;
+ }
+
+ rval = fgets(line, 1023, fh);
if ( ferror(fh) ) {
ERROR("Failed to get next filename from list.\n");
rval = NULL;
}
+
return rval;
}
-static void process_image(void *qp, void *pp, int cookie)
+static void process_image(const struct index_args *iargs,
+ struct pattern_args *pargs, int cookie)
{
- struct index_args *pargs = pp;
- struct queue_args *qargs = qp;
float *data_for_measurement;
size_t data_size;
- UnitCell *cell = qargs->static_args.cell;
- int config_cmfilter = qargs->static_args.config_cmfilter;
- int config_noisefilter = qargs->static_args.config_noisefilter;
- int config_verbose = qargs->static_args.config_verbose;
- IndexingMethod *indm = qargs->static_args.indm;
- struct beam_params *beam = qargs->static_args.beam;
+ UnitCell *cell = iargs->cell;
+ int config_cmfilter = iargs->config_cmfilter;
+ int config_noisefilter = iargs->config_noisefilter;
+ int config_verbose = iargs->config_verbose;
+ IndexingMethod *indm = iargs->indm;
+ struct beam_params *beam = iargs->beam;
int r, check;
struct hdfile *hdfile;
- char *outfile = qargs->static_args.outfile;
-
+ char *outfile = iargs->outfile;
struct image image;
+ char *outfilename;
+ int fd;
+ FILE *fh;
+ struct flock fl = {F_WRLCK, SEEK_SET, 0, 0, 0};
+
image.features = NULL;
image.data = NULL;
image.flags = NULL;
image.indexed_cell = NULL;
- image.det = copy_geom(qargs->static_args.det);
- image.copyme = qargs->static_args.copyme;
+ image.det = copy_geom(iargs->det);
+ image.copyme = iargs->copyme;
image.beam = beam;
- image.id = cookie; // MUST SET ID FOR MOSFLM TO WORK PROPERLY
+ image.id = cookie;
+ image.filename = pargs->filename;
- if (beam == NULL) {
- ERROR("Warning: no beam parameters file.\n");
- ERROR("I'm going to assume 1 ADU per photon, which is almost");
- ERROR(" certainly wrong. Peak sigmas will be incorrect.\n");
- }
+ hdfile = hdfile_open(image.filename);
+ if ( hdfile == NULL ) return;
- char *filename = NULL;
- char *imagename = pargs->filename;
- chomp(imagename);
- filename = malloc(strlen(qargs->prefix) + strlen(imagename) + 1);
- snprintf(filename, LINE_LENGTH - 1, "%s%s", qargs->prefix, imagename);
- image.filename = filename;
- hdfile = hdfile_open(filename);
- if (hdfile == NULL) return;
-
- r = hdfile_set_first_image(hdfile, "/"); // Need this to read hdf5 files
- if (r) {
+ r = hdfile_set_first_image(hdfile, "/");
+ if ( r ) {
ERROR("Couldn't select first path\n");
hdfile_close(hdfile);
return;
}
check = hdf5_read(hdfile, &image, 1);
- if (check == 1) {
+ if ( check ) {
hdfile_close(hdfile);
return;
}
- if ((image.width != image.det->max_fs + 1)
- || (image.height != image.det->max_ss + 1)) {
+ if ( (image.width != image.det->max_fs + 1 )
+ || (image.height != image.det->max_ss + 1))
+ {
ERROR("Image size doesn't match geometry size"
" - rejecting image.\n");
ERROR("Image size: %i,%i. Geometry size: %i,%i\n",
@@ -304,8 +301,8 @@ static void process_image(void *qp, void *pp, int cookie)
return;
}
- if (image.lambda < 0.0) {
- if (beam != NULL) {
+ if ( image.lambda < 0.0 ) {
+ if ( beam != NULL ) {
ERROR("Using nominal photon energy of %.2f eV\n",
beam->photon_energy);
image.lambda = ph_en_to_lambda(
@@ -320,37 +317,40 @@ static void process_image(void *qp, void *pp, int cookie)
}
fill_in_values(image.det, hdfile);
- if (config_cmfilter) {
+ if ( config_cmfilter ) {
filter_cm(&image);
}
- // Take snapshot of image after CM subtraction but before
- // the aggressive noise filter.
- data_size = image.width * image.height * sizeof (float);
+ /* Take snapshot of image after CM subtraction but before
+ * the aggressive noise filter. */
+ data_size = image.width * image.height * sizeof(float);
data_for_measurement = malloc(data_size);
- if (config_noisefilter) {
+ if ( config_noisefilter ) {
filter_noise(&image, data_for_measurement);
} else {
memcpy(data_for_measurement, image.data, data_size);
}
- switch (qargs->static_args.peaks) {
+ switch ( iargs->peaks ) {
+
case PEAK_HDF5:
// Get peaks from HDF5
- if (get_peaks(&image, hdfile,
- qargs->static_args.hdf5_peak_path)) {
- ERROR("Failed to get peaks from HDF5 file.\n");
- }
- break;
+ if (get_peaks(&image, hdfile,
+ iargs->hdf5_peak_path)) {
+ ERROR("Failed to get peaks from HDF5 file.\n");
+ }
+ break;
+
case PEAK_ZAEF:
- search_peaks(&image, qargs->static_args.threshold,
- qargs->static_args.min_gradient,
- qargs->static_args.min_snr,
- qargs->static_args.ir_inn,
- qargs->static_args.ir_mid,
- qargs->static_args.ir_out);
- break;
+ search_peaks(&image, iargs->threshold,
+ iargs->min_gradient,
+ iargs->min_snr,
+ iargs->ir_inn,
+ iargs->ir_mid,
+ iargs->ir_out);
+ break;
+
}
/* Get rid of noise-filtered version at this point
@@ -364,37 +364,32 @@ static void process_image(void *qp, void *pp, int cookie)
image.profile_radius = 0.0001e9;
/* RUN INDEXING HERE */
- index_pattern(&image, cell, indm, qargs->static_args.cellr,
- config_verbose, qargs->static_args.ipriv,
- qargs->static_args.config_insane, qargs->static_args.tols);
+ index_pattern(&image, cell, indm, iargs->cellr,
+ config_verbose, iargs->ipriv,
+ iargs->config_insane, iargs->tols);
- if (image.indexed_cell != NULL) {
+ if ( image.indexed_cell != NULL ) {
pargs->indexable = 1;
image.reflections = find_intersections(&image,
image.indexed_cell);
if (image.reflections != NULL) {
integrate_reflections(&image,
- qargs->static_args.config_closer,
- qargs->static_args.config_bgsub,
- qargs->static_args.min_int_snr,
- qargs->static_args.ir_inn,
- qargs->static_args.ir_mid,
- qargs->static_args.ir_out);
+ iargs->config_closer,
+ iargs->config_bgsub,
+ iargs->min_int_snr,
+ iargs->ir_inn,
+ iargs->ir_mid,
+ iargs->ir_out);
}
} else {
image.reflections = NULL;
}
/* Write Lock */
- struct flock fl = {F_WRLCK, SEEK_SET, 0, 0, 0};
- int fd;
fl.l_pid = getpid();
- char *outfilename = NULL;
- chomp(outfile);
- outfilename = malloc(strlen(outfile) + 1);
- snprintf(outfilename, LINE_LENGTH - 1, "%s", outfile);
- if ((fd = open(outfilename, O_WRONLY)) == -1) {
+ fd = open(outfile, O_WRONLY);
+ if ( fd == -1) {
ERROR("Error on opening\n");
exit(1);
}
@@ -404,25 +399,21 @@ static void process_image(void *qp, void *pp, int cookie)
}
/* LOCKED! Write chunk */
- FILE *fh;
fh = fopen(outfilename, "a");
if (fh == NULL) {
ERROR("Error inside lock\n");
}
- write_chunk(fh, &image, hdfile, qargs->static_args.stream_flags);
+ write_chunk(fh, &image, hdfile, iargs->stream_flags);
fclose(fh);
/* Unlock stream for other processes */
fl.l_type = F_UNLCK; /* set to unlock same region */
- if (fcntl(fd, F_SETLK, &fl) == -1) {
+ if ( fcntl(fd, F_SETLK, &fl) == -1 ) {
ERROR("fcntl");
exit(1);
}
close(fd);
- qargs->n_indexable += pargs->indexable;
- qargs->n_processed++;
-
/* Only free cell if found */
cell_free(image.indexed_cell);
@@ -436,25 +427,42 @@ static void process_image(void *qp, void *pp, int cookie)
static void run_work(const struct index_args *iargs,
- int filename_pipe, int results_pipe)
+ int filename_pipe, int results_pipe, int cookie)
{
int allDone = 0;
while ( !allDone ) {
- /* read from pipe and return number of bytes read */
- if ((buff_count=read(fd_pipeOut[batchNum-1][0],&buffR,BUFFER))<0) {
- ERROR("read1");
- } else if (buff_count > 0) {
- /* process image */
- pargs.filename = buffR;
+ struct pattern_args pargs;
+ int r, w;
+ char buf[1024];
+
+ r = read(filename_pipe, buf, 1024);
+
+ if ( r < 0 ) {
+
+ ERROR("read() failed!\n");
+
+ } else if ( r > 0 ) {
+
+ int c;
+
+ /* Process image */
+ pargs.filename = buf;
pargs.indexable = 0;
- process_image(&qargs, &pargs, batchNum);
- /* request another image */
- buff_count = sprintf(buffW, "%d\n", pargs.indexable);
- if(write (fd_pipeIn[batchNum-1][1], buffW, buff_count)<0)
+
+ STATUS("Got filename: '%s'\n", buf);
+
+ process_image(iargs, &pargs, cookie);
+
+ /* Request another image */
+ c = sprintf(buf, "%i\n", pargs.indexable);
+ w = write(results_pipe, buf, c);
+ if ( w < 0 ) {
ERROR("write P0");
- } else if (buff_count == 0) {
+ }
+
+ } else {
allDone = 1;
}
@@ -553,14 +561,12 @@ int main(int argc, char *argv[])
int peaks;
int n_proc = 1;
char *prepare_line;
- char prepare_filename[LINE_LENGTH];
- struct queue_args qargs;
- struct index_args pargs;
+ char prepare_filename[1024];
+ struct index_args iargs;
struct beam_params *beam = NULL;
char *element = NULL;
double nominal_photon_energy;
int stream_flags = STREAM_INTEGRATED;
- char *endptr;
char *hdf5_peak_path = NULL;
struct copy_hdf5_field *copyme;
char *intrad = NULL;
@@ -573,6 +579,9 @@ int main(int argc, char *argv[])
pid_t *pids;
int *filename_pipes;
int *result_pipes;
+ fd_set fds;
+ int i;
+ int allDone, nFinished;
copyme = new_copy_hdf5_field_list();
if ( copyme == NULL ) {
@@ -793,7 +802,7 @@ int main(int argc, char *argv[])
}
}
- if ( nProcesses == 0 ) {
+ if ( n_proc == 0 ) {
ERROR("Invalid number of processes.\n");
return 1;
}
@@ -900,18 +909,15 @@ int main(int argc, char *argv[])
} else {
STATUS("No beam parameters file was given, so I'm taking the"
" nominal photon energy to be 2 keV.\n");
+ ERROR("I'm also going to assume 1 ADU per photon, which is");
+ ERROR(" almost certainly wrong. Peak sigmas will be"
+ " incorrect.\n");
nominal_photon_energy = 2000.0;
}
- if ( beam == NULL ) {
- ERROR("Warning: no beam parameters file.\n");
- ERROR("I'm going to assume 1 ADU per photon, which is almost");
- ERROR(" certainly wrong. Peak sigmas will be incorrect.\n");
- }
-
/* Get first filename and use it to set up the indexing */
- prepare_line = malloc(LINE_LENGTH*sizeof(char));
- rval = fgets(prepare_line, LINE_LENGTH-1, fh);
+ prepare_line = malloc(1024);
+ rval = fgets(prepare_line, 1023, fh);
if ( rval == NULL ) {
ERROR("Failed to get filename to prepare indexing.\n");
return 1;
@@ -923,7 +929,7 @@ int main(int argc, char *argv[])
free(prepare_line);
prepare_line = tmp;
}
- snprintf(prepare_filename, LINE_LENGTH-1, "%s%s", prefix, prepare_line);
+ snprintf(prepare_filename, 1023, "%s%s", prefix, prepare_line);
rewind(fh);
/* Prepare the indexer */
@@ -979,6 +985,20 @@ int main(int argc, char *argv[])
n_processed_last_stats = 0;
t_last_stats = get_monotonic_seconds();
+ FD_ZERO(&fds);
+ filename_pipes = calloc(n_proc, sizeof(int));
+ result_pipes = calloc(n_proc, sizeof(int));
+ if ( (filename_pipes == NULL) || (result_pipes == NULL) ) {
+ ERROR("Couldn't allocate memory for pipes.\n");
+ return 1;
+ }
+
+ pids = calloc(n_proc, sizeof(pid_t));
+ if ( pids == NULL ) {
+ ERROR("Couldn't allocate memory for PIDs.\n");
+ return 1;
+ }
+
/* Fork the right number of times */
for ( i=0; i<n_proc; i++ ) {
@@ -991,7 +1011,7 @@ int main(int argc, char *argv[])
return 1;
}
- if ( pipe(results_pipe) == - 1 ) {
+ if ( pipe(result_pipe) == - 1 ) {
ERROR("pipe() failed!\n");
return 1;
}
@@ -1007,80 +1027,116 @@ int main(int argc, char *argv[])
* pipe, and the 'write' end of the result pipe. */
close(filename_pipe[1]);
close(result_pipe[0]);
- run_work(&iargs, filename_pipe[0], result_pipe[1]);
+ run_work(&iargs, filename_pipe[0], result_pipe[1], i);
exit(0);
}
/* Parent process gets the 'write' end of the filename pipe
* and the 'read' end of the result pipe. */
- pid[i] = p;
+ pids[i] = p;
close(filename_pipe[0]);
close(result_pipe[1]);
filename_pipes[i] = filename_pipe[1];
result_pipes[i] = result_pipe[0];
+ FD_SET(result_pipes[i], &fds);
+
}
/* Send first image to all children */
- char *nextImage = NULL;
- for ( i=0; i<nProcesses; i++ ) {
+ for ( i=0; i<n_proc; i++ ) {
+
+ char *nextImage;
+
nextImage = get_pattern(fh);
- buff_count = sprintf(buffW, "%s",nextImage);
- write (fd_pipeOut[i][1], buffW, buff_count);
+
+ write(filename_pipes[i], nextImage, strlen(nextImage));
+ write(filename_pipes[i], "\n", 1);
+
}
- int nFinished = 0;
- while (!allDone) {
- /* select from file set for reading */
- if ((ready_fd = select(max_fd,&fdset,NULL,NULL,NULL)) < 0)
- ERROR("select");
- if (ready_fd > 0) {
- for ( i=0; i<nProcesses; i++ ) {
- /* is in file set that raised flag? */
- if (FD_ISSET(fd_pipeIn[i][0],&fdset)) {
- /* read from pipe and return number of bytes read */
- if ((buff_count=read(fd_pipeIn[i][0],&buffR,BUFFER))<0) {
- ERROR("read");
- } else {
- qargs.n_indexable += atoi(buffR);
- qargs.n_processed++;
- /* write to pipe */
- if ((nextImage = get_pattern(fh)) == NULL){
- nFinished++; /* no more images */
- if ( nFinished == nProcesses )
- allDone = 1; /* EXIT */
- } else {
- /* send out image */
- buff_count = sprintf(buffW, "%s",nextImage);
- if (write (fd_pipeOut[i][1], buffW, buff_count)<0)
- ERROR("write pipe");
- }
+
+ allDone = 0;
+ nFinished = 0;
+ while ( !allDone ) {
+
+ int r;
+ struct timeval tv;
+ fd_set fds_copy;
+ double tNow;
+
+ tv.tv_sec = 5;
+ tv.tv_usec = 0;
+
+ memcpy(&fds_copy, &fds, sizeof(fd_set));
+ r = select(n_proc, &fds, NULL, NULL, &tv);
+
+ if ( r < 0 ) {
+
+ ERROR("select() failed!\n");
+
+ } else {
+
+ for ( i=0; i<n_proc; i++ ) {
+
+ char *nextImage;
+ char results[1024];
+
+ if ( !FD_ISSET(result_pipes[i], &fds_copy) ) {
+ continue;
+ }
+
+ r = read(result_pipes[i], results, 1024);
+ if ( r < 0 ) {
+ ERROR("read() failed!");
+ continue;
+ }
+
+ n_indexable += atoi(results);
+ n_processed++;
+
+ /* Send next filename */
+ nextImage = get_pattern(fh);
+ if ( nextImage == NULL ) {
+ /* no more images */
+ nFinished++;
+ if ( nFinished == n_proc ) {
+ allDone = 1;
+ }
+ } else {
+
+ r = write(filename_pipes[i], nextImage,
+ strlen(nextImage));
+ if ( r < 0 ) {
+ ERROR("write pipe");
}
}
+
}
+
}
- /* file set is modified, so copy original from tmpset */
- memcpy((void *) &fdset,(void *) &tmpset, sizeof(fd_set));
-
- /* Update to screen */
- double tNow = get_monotonic_seconds();
- if ( tNow >= qargs.t_last_stats+STATS_EVERY_N_SECONDS ) {
- STATUS("%i out of %i indexed so far,"
- " %i out of %i since the last message.\n\n",
- qargs.n_indexable, qargs.n_processed,
- qargs.n_indexable - qargs.n_indexable_last_stats,
- qargs.n_processed - qargs.n_processed_last_stats);
-
- qargs.n_indexable_last_stats = qargs.n_indexable;
- qargs.n_processed_last_stats = qargs.n_processed;
- qargs.t_last_stats = tNow;
+
+ /* Update progress */
+ tNow = get_monotonic_seconds();
+ if ( tNow >= t_last_stats+STATS_EVERY_N_SECONDS ) {
+
+ STATUS("%i out of %i indexed so far,"
+ " %i out of %i since the last message.\n\n",
+ n_indexable, n_processed,
+ n_indexable - n_indexable_last_stats,
+ n_processed - n_processed_last_stats);
+
+ n_indexable_last_stats = n_indexable;
+ n_processed_last_stats = n_processed;
+ t_last_stats = tNow;
+
}
+
}
- /* close my pipes */
- for ( i=0; i<nProcesses; i++ ) {
- close(fd_pipeIn[i][0]);
- close(fd_pipeOut[i][1]);
+
+ for ( i=0; i<n_proc; i++ ) {
+ close(filename_pipes[i]);
+ close(result_pipes[i]);
}
- tEnd = get_monotonic_seconds();
cleanup_indexing(ipriv);
@@ -1096,9 +1152,8 @@ int main(int argc, char *argv[])
if ( fh != stdin ) fclose(fh);
if ( ofh != stdout ) fclose(ofh);
- if (batchNum == qargs.updateReader) {
- STATUS("There were %i images, of which %i could be indexed.\n",
- qargs.n_processed, qargs.n_indexable);
- }
+ STATUS("There were %i images, of which %i could be indexed.\n",
+ n_processed, n_indexable);
+
return 0;
}