aboutsummaryrefslogtreecommitdiff
path: root/libcrystfel
diff options
context:
space:
mode:
Diffstat (limited to 'libcrystfel')
-rw-r--r--libcrystfel/src/crystal.c11
-rw-r--r--libcrystfel/src/crystal.h9
-rw-r--r--libcrystfel/src/image.h2
-rw-r--r--libcrystfel/src/stream.c451
-rw-r--r--libcrystfel/src/stream.h32
5 files changed, 294 insertions, 211 deletions
diff --git a/libcrystfel/src/crystal.c b/libcrystfel/src/crystal.c
index 85430b6f..3c0649cc 100644
--- a/libcrystfel/src/crystal.c
+++ b/libcrystfel/src/crystal.c
@@ -50,10 +50,8 @@
struct _crystal
{
/* Information about the crystal */
- UnitCell *indexed_cell;
- UnitCell *candidate_cells[MAX_CELL_CANDIDATES];
- int ncells;
- double m; /* Mosaicity in radians */
+ UnitCell *cell;
+ double m; /* Mosaicity in radians */
double osf;
double profile_radius;
int pr_dud;
@@ -83,6 +81,11 @@ Crystal *crystal_new()
cryst = malloc(sizeof(Crystal));
if ( cryst == NULL ) return NULL;
+ cryst->cell = NULL;
+ cryst->reflections = NULL;
+ cryst->diffracting_resolution = 0.0;
+ cryst->n_saturated = 0;
+
return cryst;
}
diff --git a/libcrystfel/src/crystal.h b/libcrystfel/src/crystal.h
index 466116a4..e5e40589 100644
--- a/libcrystfel/src/crystal.h
+++ b/libcrystfel/src/crystal.h
@@ -53,5 +53,14 @@ extern void crystal_free(Crystal *cryst);
extern UnitCell *crystal_get_cell(Crystal *cryst);
extern double crystal_get_profile_radius(Crystal *cryst);
extern RefList *crystal_get_reflections(Crystal *cryst);
+extern double crystal_get_resolution_limit(Crystal *cryst);
+extern long long int crystal_get_num_saturated_reflections(Crystal *cryst);
+
+extern void crystal_set_cell(Crystal *cryst, UnitCell *cell);
+extern void crystal_set_profile_radius(Crystal *cryst, double r);
+extern void crystal_set_reflections(Crystal *cryst, RefList *reflist);
+extern void crystal_set_resolution_limit(Crystal *cryst, double res);
+extern void crystal_set_num_saturated_reflections(Crystal *cryst, int n);
+
#endif /* CRYSTAL_H */
diff --git a/libcrystfel/src/image.h b/libcrystfel/src/image.h
index bd90a920..d54cf43f 100644
--- a/libcrystfel/src/image.h
+++ b/libcrystfel/src/image.h
@@ -132,7 +132,7 @@ struct image {
uint16_t *flags;
double *twotheta;
- struct crystal **crystals;
+ Crystal **crystals;
int n_crystals;
struct detector *det;
diff --git a/libcrystfel/src/stream.c b/libcrystfel/src/stream.c
index b8ff9238..943ee656 100644
--- a/libcrystfel/src/stream.c
+++ b/libcrystfel/src/stream.c
@@ -3,12 +3,12 @@
*
* Stream tools
*
- * Copyright © 2012 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2013 Deutsches Elektronen-Synchrotron DESY,
* a research centre of the Helmholtz Association.
* Copyright © 2012 Richard Kirian
*
* Authors:
- * 2010-2012 Thomas White <taw@physics.org>
+ * 2010-2013 Thomas White <taw@physics.org>
* 2011 Richard Kirian
* 2011 Andrew Aquila
*
@@ -37,6 +37,7 @@
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
+#include <assert.h>
#include "cell.h"
#include "utils.h"
@@ -45,83 +46,27 @@
#include "reflist.h"
#include "reflist-utils.h"
+#define LATEST_MAJOR_VERSION (2)
+#define LATEST_MINOR_VERSION (1)
#define CHUNK_START_MARKER "----- Begin chunk -----"
#define CHUNK_END_MARKER "----- End chunk -----"
#define PEAK_LIST_START_MARKER "Peaks from peak search"
#define PEAK_LIST_END_MARKER "End of peak list"
+#define CRYSTAL_START_MARKER "--- Begin crystal"
+#define CRYSTAL_END_MARKER "--- End crystal"
#define REFLECTION_START_MARKER "Reflections measured after indexing"
/* REFLECTION_END_MARKER is over in reflist-utils.h because it is also
* used to terminate a standalone list of reflections */
-static void exclusive(const char *a, const char *b)
-{
- ERROR("The stream options '%s' and '%s' are mutually exclusive.\n",
- a, b);
-}
-
-int parse_stream_flags(const char *a)
+struct _stream
{
- int n, i;
- int ret = STREAM_NONE;
- char **flags;
-
- n = assplode(a, ",", &flags, ASSPLODE_NONE);
-
- for ( i=0; i<n; i++ ) {
-
- if ( strcmp(flags[i], "integrated") == 0) {
-
- ret |= STREAM_INTEGRATED;
-
- } else if ( strcmp(flags[i], "peaks") == 0) {
- if ( ret & STREAM_PEAKS_IF_INDEXED ) {
- exclusive("peaks", "peaksifindexed");
- return -1;
- }
- if ( ret & STREAM_PEAKS_IF_NOT_INDEXED ) {
- exclusive("peaks", "peaksifnotindexed");
- return -1;
- }
- ret |= STREAM_PEAKS;
-
- } else if ( strcmp(flags[i], "peaksifindexed") == 0) {
- if ( ret & STREAM_PEAKS ) {
- exclusive("peaks", "peaksifindexed");
- return -1;
- }
- if ( ret & STREAM_PEAKS_IF_NOT_INDEXED ) {
- exclusive("peaksifnotindexed",
- "peaksifindexed");
- return -1;
- }
- ret |= STREAM_PEAKS_IF_INDEXED;
-
- } else if ( strcmp(flags[i], "peaksifnotindexed") == 0) {
- if ( ret & STREAM_PEAKS ) {
- exclusive("peaks", "peaksifnotindexed");
- return -1;
- }
- if ( ret & STREAM_PEAKS_IF_INDEXED ) {
- exclusive("peaksifnotindexed",
- "peaksifindexed");
- return -1;
- }
- ret |= STREAM_PEAKS_IF_NOT_INDEXED;
-
- } else {
- ERROR("Unrecognised stream flag '%s'\n", flags[i]);
- return -1;
- }
-
- free(flags[i]);
-
- }
- free(flags);
+ FILE *fh;
- return ret;
-}
+ int major_version;
+ int minor_version;
+};
int count_patterns(FILE *fh)
@@ -215,96 +160,106 @@ static void write_peaks(struct image *image, FILE *ofh)
}
-void write_chunk(FILE *ofh, struct image *i, struct hdfile *hdfile, int f)
+static void write_crystal(Stream *st, Crystal *cr, int include_reflections)
{
+ UnitCell *cell;
+ RefList *reflist;
double asx, asy, asz;
double bsx, bsy, bsz;
double csx, csy, csz;
double a, b, c, al, be, ga;
- fprintf(ofh, CHUNK_START_MARKER"\n");
+ fprintf(st->fh, CRYSTAL_START_MARKER"\n");
- fprintf(ofh, "Image filename: %s\n", i->filename);
+ cell = crystal_get_cell(cr);
+ assert(cell != NULL);
- if ( i->indexed_cell != NULL ) {
+ cell_get_parameters(cell, &a, &b, &c, &al, &be, &ga);
+ fprintf(st->fh, "Cell parameters %7.5f %7.5f %7.5f nm,"
+ " %7.5f %7.5f %7.5f deg\n",
+ a*1.0e9, b*1.0e9, c*1.0e9,
+ rad2deg(al), rad2deg(be), rad2deg(ga));
- cell_get_parameters(i->indexed_cell, &a, &b, &c,
- &al, &be, &ga);
- fprintf(ofh, "Cell parameters %7.5f %7.5f %7.5f nm,"
- " %7.5f %7.5f %7.5f deg\n",
- a*1.0e9, b*1.0e9, c*1.0e9,
- rad2deg(al), rad2deg(be), rad2deg(ga));
+ cell_get_reciprocal(cell, &asx, &asy, &asz,
+ &bsx, &bsy, &bsz,
+ &csx, &csy, &csz);
+ fprintf(st->fh, "astar = %+9.7f %+9.7f %+9.7f nm^-1\n",
+ asx/1e9, asy/1e9, asz/1e9);
+ fprintf(st->fh, "bstar = %+9.7f %+9.7f %+9.7f nm^-1\n",
+ bsx/1e9, bsy/1e9, bsz/1e9);
+ fprintf(st->fh, "cstar = %+9.7f %+9.7f %+9.7f nm^-1\n",
+ csx/1e9, csy/1e9, csz/1e9);
- cell_get_reciprocal(i->indexed_cell, &asx, &asy, &asz,
- &bsx, &bsy, &bsz,
- &csx, &csy, &csz);
- fprintf(ofh, "astar = %+9.7f %+9.7f %+9.7f nm^-1\n",
- asx/1e9, asy/1e9, asz/1e9);
- fprintf(ofh, "bstar = %+9.7f %+9.7f %+9.7f nm^-1\n",
- bsx/1e9, bsy/1e9, bsz/1e9);
- fprintf(ofh, "cstar = %+9.7f %+9.7f %+9.7f nm^-1\n",
- csx/1e9, csy/1e9, csz/1e9);
+ reflist = crystal_get_reflections(cr);
+ if ( reflist != NULL ) {
- } else {
+ fprintf(st->fh, "diffraction_resolution_limit"
+ " = %.2f nm^-1 or %.2f A\n",
+ crystal_get_resolution_limit(cr)/1e9,
+ 1e9 / crystal_get_resolution_limit(cr));
- fprintf(ofh, "No unit cell from indexing.\n");
+ fprintf(st->fh, "num_saturated_reflections = %lli\n",
+ crystal_get_num_saturated_reflections(cr));
}
- fprintf(ofh, "photon_energy_eV = %f\n",
- J_to_eV(ph_lambda_to_en(i->lambda)));
+ if ( include_reflections ) {
+
+ fprintf(st->fh, "\n");
- if ( i->reflections != NULL ) {
+ if ( reflist != NULL ) {
- fprintf(ofh, "diffraction_resolution_limit"
- " = %.2f nm^-1 or %.2f A\n",
- i->diffracting_resolution/1e9,
- 1e9 / i->diffracting_resolution);
+ fprintf(st->fh, REFLECTION_START_MARKER"\n");
+ write_reflections_to_file(st->fh, reflist);
+ fprintf(st->fh, REFLECTION_END_MARKER"\n");
+
+ } else {
- fprintf(ofh, "num_saturated_reflections"
- " = %lli\n", i->n_saturated);
+ fprintf(st->fh, "No integrated reflections.\n");
+ }
}
+ fprintf(st->fh, CRYSTAL_START_MARKER"\n\n");
+}
+
+
+void write_chunk(Stream *st, struct image *i, struct hdfile *hdfile,
+ int include_peaks, int include_reflections)
+{
+ int j;
+
+ fprintf(st->fh, CHUNK_START_MARKER"\n");
+
+ fprintf(st->fh, "Image filename: %s\n", i->filename);
+
if ( i->det != NULL ) {
int j;
for ( j=0; j<i->det->n_panels; j++ ) {
- fprintf(ofh, "camera_length_%s = %f\n",
+ fprintf(st->fh, "camera_length_%s = %f\n",
i->det->panels[j].name, i->det->panels[j].clen);
}
}
- copy_hdf5_fields(hdfile, i->copyme, ofh);
+ copy_hdf5_fields(hdfile, i->copyme, st->fh);
- if ( (f & STREAM_PEAKS)
- || ((f & STREAM_PEAKS_IF_INDEXED) && (i->indexed_cell != NULL))
- || ((f & STREAM_PEAKS_IF_NOT_INDEXED) && (i->indexed_cell == NULL)) )
- {
- fprintf(ofh, "\n");
- write_peaks(i, ofh);
+ if ( include_peaks ) {
+ fprintf(st->fh, "\n");
+ write_peaks(i, st->fh);
}
- if ( f & STREAM_INTEGRATED ) {
-
- fprintf(ofh, "\n");
-
- if ( i->reflections != NULL ) {
-
- fprintf(ofh, REFLECTION_START_MARKER"\n");
- write_reflections_to_file(ofh, i->reflections);
- fprintf(ofh, REFLECTION_END_MARKER"\n");
-
- } else {
-
- fprintf(ofh, "No integrated reflections.\n");
+ fprintf(st->fh, "photon_energy_eV = %f\n",
+ J_to_eV(ph_lambda_to_en(i->lambda)));
- }
+ fprintf(st->fh, "\n");
+ for ( j=0; j<i->n_crystals; j++ ) {
+ write_crystal(st, i->crystals[j], include_reflections);
}
- fprintf(ofh, CHUNK_END_MARKER"\n\n");
+ fprintf(st->fh, CHUNK_END_MARKER"\n\n");
}
@@ -328,8 +283,7 @@ static int find_start_of_chunk(FILE *fh)
}
-/* Read the next chunk from a stream and fill in 'image' */
-int read_chunk(FILE *fh, struct image *image)
+void read_crystal(Stream *st, struct image *image)
{
char line[1024];
char *rval = NULL;
@@ -337,22 +291,116 @@ int read_chunk(FILE *fh, struct image *image)
int have_as = 0;
int have_bs = 0;
int have_cs = 0;
+ Crystal *cr;
+ int n;
+ Crystal **crystals_new;
+
+ cr = crystal_new();
+ if ( cr == NULL ) {
+ ERROR("Failed to allocate crystal!\n");
+ return;
+ }
+
+ do {
+
+ float u, v, w;
+
+ rval = fgets(line, 1023, st->fh);
+
+ /* Trouble? */
+ if ( rval == NULL ) break;
+
+ chomp(line);
+ if ( sscanf(line, "astar = %f %f %f", &u, &v, &w) == 3 ) {
+ as.u = u*1e9; as.v = v*1e9; as.w = w*1e9;
+ have_as = 1;
+ }
+
+ if ( sscanf(line, "bstar = %f %f %f", &u, &v, &w) == 3 ) {
+ bs.u = u*1e9; bs.v = v*1e9; bs.w = w*1e9;
+ have_bs = 1;
+ }
+
+ if ( sscanf(line, "cstar = %f %f %f", &u, &v, &w) == 3 ) {
+ cs.u = u*1e9; cs.v = v*1e9; cs.w = w*1e9;
+ have_cs = 1;
+ }
+
+ if ( have_as && have_bs && have_cs ) {
+
+ UnitCell *cell;
+
+ cell = crystal_get_cell(cr);
+
+ if ( cell != NULL ) {
+ ERROR("Duplicate cell found in stream!\n");
+ ERROR("I'll use the most recent one.\n");
+ cell_free(cell);
+ }
+
+ cell = cell_new_from_reciprocal_axes(as, bs, cs);
+ crystal_set_cell(cr, cell);
+
+ have_as = 0; have_bs = 0; have_cs = 0;
+
+ }
+
+ if ( strncmp(line, "num_saturated_reflections = ", 28) == 0 ) {
+ int n = atoi(line+28);
+ crystal_set_num_saturated_reflections(cr, n);
+ }
+
+ if ( strcmp(line, REFLECTION_START_MARKER) == 0 ) {
+
+ RefList *reflections;
+
+ reflections = read_reflections_from_file(st->fh);
+
+ if ( reflections == NULL ) {
+ ERROR("Failed while reading reflections\n");
+ break;
+ }
+
+ crystal_set_reflections(cr, reflections);
+
+ }
+
+ if ( strcmp(line, CRYSTAL_END_MARKER) == 0 ) break;
+
+ } while ( 1 );
+
+ /* Add crystal to the list for this image */
+ n = image->n_crystals+1;
+ crystals_new = realloc(image->crystals, n*sizeof(Crystal *));
+
+ if ( crystals_new == NULL ) {
+ ERROR("Failed to expand crystal list!\n");
+ } else {
+ image->crystals = crystals_new;
+ image->crystals[image->n_crystals++] = cr;
+ }
+
+}
+
+
+/* Read the next chunk from a stream and fill in 'image' */
+int read_chunk(Stream *st, struct image *image)
+{
+ char line[1024];
+ char *rval = NULL;
int have_filename = 0;
int have_ev = 0;
- if ( find_start_of_chunk(fh) ) return 1;
+ if ( find_start_of_chunk(st->fh) ) return 1;
image->lambda = -1.0;
image->features = NULL;
- image->reflections = NULL;
- image->indexed_cell = NULL;
- image->n_saturated = 0;
+ image->crystals = NULL;
+ image->n_crystals = 0;
do {
- float u, v, w;
-
- rval = fgets(line, 1023, fh);
+ rval = fgets(line, 1023, st->fh);
/* Trouble? */
if ( rval == NULL ) break;
@@ -364,6 +412,11 @@ int read_chunk(FILE *fh, struct image *image)
have_filename = 1;
}
+ if ( strncmp(line, "photon_energy_eV = ", 19) == 0 ) {
+ image->lambda = ph_en_to_lambda(eV_to_J(atof(line+19)));
+ have_ev = 1;
+ }
+
if ( strncmp(line, "camera_length_", 14) == 0 ) {
if ( image->det != NULL ) {
@@ -390,56 +443,19 @@ int read_chunk(FILE *fh, struct image *image)
}
}
- if ( sscanf(line, "astar = %f %f %f", &u, &v, &w) == 3 ) {
- as.u = u*1e9; as.v = v*1e9; as.w = w*1e9;
- have_as = 1;
- }
-
- if ( sscanf(line, "bstar = %f %f %f", &u, &v, &w) == 3 ) {
- bs.u = u*1e9; bs.v = v*1e9; bs.w = w*1e9;
- have_bs = 1;
- }
-
- if ( sscanf(line, "cstar = %f %f %f", &u, &v, &w) == 3 ) {
- cs.u = u*1e9; cs.v = v*1e9; cs.w = w*1e9;
- have_cs = 1;
- }
-
- if ( have_as && have_bs && have_cs ) {
- if ( image->indexed_cell != NULL ) {
- ERROR("Duplicate cell found in stream!\n");
- cell_free(image->indexed_cell);
- }
- image->indexed_cell = cell_new_from_reciprocal_axes(as,
- bs,
- cs);
- have_as = 0; have_bs = 0; have_cs = 0;
- }
-
- if ( strncmp(line, "photon_energy_eV = ", 19) == 0 ) {
- image->lambda = ph_en_to_lambda(eV_to_J(atof(line+19)));
- have_ev = 1;
- }
-
- if ( strncmp(line, "num_saturated_reflections = ", 28) == 0 ) {
- image->n_saturated = atoi(line+28);
- }
-
if ( strcmp(line, PEAK_LIST_START_MARKER) == 0 ) {
- if ( read_peaks(fh, image) ) {
+ if ( read_peaks(st->fh, image) ) {
ERROR("Failed while reading peaks\n");
return 1;
}
}
- if ( strcmp(line, REFLECTION_START_MARKER) == 0 ) {
- image->reflections = read_reflections_from_file(fh);
- if ( image->reflections == NULL ) {
- ERROR("Failed while reading reflections\n");
- return 1;
- }
+ if ( strcmp(line, CRYSTAL_START_MARKER) == 0 ) {
+ read_crystal(st, image);
}
+ /* A chunk must have at least a filename and a wavelength,
+ * otherwise it's incomplete */
if ( strcmp(line, CHUNK_END_MARKER) == 0 ) {
if ( have_filename && have_ev ) return 0;
ERROR("Incomplete chunk found in input file.\n");
@@ -457,7 +473,6 @@ void write_stream_header(FILE *ofh, int argc, char *argv[])
{
int i;
- fprintf(ofh, "CrystFEL stream format 2.0\n");
fprintf(ofh, "Command line:");
for ( i=0; i<argc; i++ ) {
fprintf(ofh, " %s", argv[i]);
@@ -489,24 +504,90 @@ int skip_some_files(FILE *fh, int n)
}
-int is_stream(const char *filename)
+Stream *open_stream_for_read(const char *filename)
{
- FILE *fh;
- char line[1024];
- char *rval = NULL;
- fh = fopen(filename, "r");
- if ( fh == NULL ) {
- return -1;
+ Stream *st;
+
+ st = malloc(sizeof(struct _stream));
+ if ( st == NULL ) return NULL;
+
+ st->fh = fopen(filename, "r");
+ if ( st->fh == NULL ) {
+ free(st);
+ return NULL;
}
- rval = fgets(line, 1023, fh);
- fclose(fh);
+
+ char line[1024];
+ char *rval;
+
+ rval = fgets(line, 1023, st->fh);
if ( rval == NULL ) {
- return -1;
+ ERROR("Failed to read stream version.\n");
+ close_stream(st);
+ return NULL;
}
+
if ( strncmp(line, "CrystFEL stream format 2.0", 26) == 0 ) {
- return 1;
+ st->major_version = 2;
+ st->minor_version = 0;
+ } else if ( strncmp(line, "CrystFEL stream format 2.1", 26) == 0 ) {
+ st->major_version = 2;
+ st->minor_version = 1;
+ } else {
+ ERROR("Invalid stream, or stream format is too new.\n");
+ close_stream(st);
+ return NULL;
}
- else {
- return 0;
+
+ return st;
+}
+
+
+Stream *open_stream_for_write(const char *filename)
+{
+ Stream *st;
+
+ st = malloc(sizeof(struct _stream));
+ if ( st == NULL ) return NULL;
+
+ st->fh = fopen(filename, "w");
+ if ( st->fh == NULL ) {
+ free(st);
+ return NULL;
}
+
+ st->major_version = LATEST_MAJOR_VERSION;
+ st->minor_version = LATEST_MINOR_VERSION;
+
+ fprintf(st->fh, "CrystFEL stream format %i.%i\n",
+ st->major_version, st->minor_version);
+
+ return st;
+}
+
+
+void close_stream(Stream *st)
+{
+ fclose(st->fh);
+ free(st);
+}
+
+
+int is_stream(const char *filename)
+{
+ FILE *fh;
+ char line[1024];
+ char *rval;
+
+ fh = fopen(filename, "r");
+ if ( fh == NULL ) return 0;
+
+ rval = fgets(line, 1023, fh);
+ fclose(fh);
+ if ( rval == NULL ) return 0;
+
+ if ( strncmp(line, "CrystFEL stream format 2.0", 26) == 0 ) return 1;
+ if ( strncmp(line, "CrystFEL stream format 2.1", 26) == 0 ) return 1;
+
+ return 0;
}
diff --git a/libcrystfel/src/stream.h b/libcrystfel/src/stream.h
index a544798f..db74583d 100644
--- a/libcrystfel/src/stream.h
+++ b/libcrystfel/src/stream.h
@@ -3,11 +3,11 @@
*
* Stream tools
*
- * Copyright © 2012 Deutsches Elektronen-Synchrotron DESY,
+ * Copyright © 2013 Deutsches Elektronen-Synchrotron DESY,
* a research centre of the Helmholtz Association.
*
* Authors:
- * 2010-2012 Thomas White <taw@physics.org>
+ * 2010-2013 Thomas White <taw@physics.org>
* 2011 Andrew Aquila
*
* This file is part of CrystFEL.
@@ -38,30 +38,20 @@
struct image;
struct hdfile;
-/* Possible options dictating what goes into the output stream */
-enum
-{
- STREAM_NONE = 0,
- STREAM_INTEGRATED = 1<<0,
- STREAM_PEAKS = 1<<2,
- STREAM_PEAKS_IF_INDEXED = 1<<3,
- STREAM_PEAKS_IF_NOT_INDEXED = 1<<4,
-};
+typedef struct _stream Stream;
+extern Stream *open_stream_for_read(const char *filename);
+extern Stream *open_stream_for_write(const char *filename);
+extern void close_stream(Stream *st);
-extern int count_patterns(FILE *fh);
-
-extern void write_stream_header(FILE *ofh, int argc, char *argv[]);
-
-extern void write_chunk(FILE *ofh, struct image *image, struct hdfile *hdfile,
- int flags);
+extern void write_chunk(Stream *st, struct image *image, struct hdfile *hdfile,
+ int include_peaks, int include_reflections);
-extern int parse_stream_flags(const char *a);
-
-extern int read_chunk(FILE *fh, struct image *image);
+extern int read_chunk(Stream *st, struct image *image);
+/* Nasty functions that should be avoided */
+extern int count_patterns(FILE *fh);
extern int skip_some_files(FILE *fh, int n);
-
extern int is_stream(const char *filename);
#endif /* STREAM_H */
generated by cgit v1.2.3 (git 2.25.1) at 2024-11-16 19:27:07 +0000