diff options
-rw-r--r-- | src/Makefile.am | 2 | ||||
-rw-r--r-- | src/Makefile.in | 9 | ||||
-rw-r--r-- | src/cell.c | 42 | ||||
-rw-r--r-- | src/cell.h | 2 | ||||
-rw-r--r-- | src/diffraction-gpu.c | 6 | ||||
-rw-r--r-- | src/diffraction-gpu.h | 4 | ||||
-rw-r--r-- | src/diffraction.c | 8 | ||||
-rw-r--r-- | src/diffraction.h | 2 | ||||
-rw-r--r-- | src/get_hkl.c | 2 | ||||
-rw-r--r-- | src/image.h | 1 | ||||
-rw-r--r-- | src/index.c | 7 | ||||
-rw-r--r-- | src/index.h | 7 | ||||
-rw-r--r-- | src/indexamajig.c | 53 | ||||
-rw-r--r-- | src/pattern_sim.c | 19 |
14 files changed, 104 insertions, 60 deletions
diff --git a/src/Makefile.am b/src/Makefile.am index e696fcbf..ef9338c5 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -23,7 +23,7 @@ process_hkl_LDADD = @LIBS@ indexamajig_SOURCES = indexamajig.c hdf5-file.c utils.c cell.c image.c \ intensities.c peaks.c index.c filters.c \ - diffraction.c detector.c sfac.c dirax.c + diffraction.c detector.c sfac.c dirax.c reflections.c indexamajig_LDADD = @LIBS@ if HAVE_OPENCL indexamajig_SOURCES += diffraction-gpu.c cl-utils.c diff --git a/src/Makefile.in b/src/Makefile.in index 87fb7e4c..4b59a13d 100644 --- a/src/Makefile.in +++ b/src/Makefile.in @@ -69,15 +69,16 @@ hdfsee_OBJECTS = $(am_hdfsee_OBJECTS) hdfsee_DEPENDENCIES = am__indexamajig_SOURCES_DIST = indexamajig.c hdf5-file.c utils.c \ cell.c image.c intensities.c peaks.c index.c filters.c \ - diffraction.c detector.c sfac.c dirax.c diffraction-gpu.c \ - cl-utils.c + diffraction.c detector.c sfac.c dirax.c reflections.c \ + diffraction-gpu.c cl-utils.c @HAVE_OPENCL_TRUE@am__objects_1 = diffraction-gpu.$(OBJEXT) \ @HAVE_OPENCL_TRUE@ cl-utils.$(OBJEXT) am_indexamajig_OBJECTS = indexamajig.$(OBJEXT) hdf5-file.$(OBJEXT) \ utils.$(OBJEXT) cell.$(OBJEXT) image.$(OBJEXT) \ intensities.$(OBJEXT) peaks.$(OBJEXT) index.$(OBJEXT) \ filters.$(OBJEXT) diffraction.$(OBJEXT) detector.$(OBJEXT) \ - sfac.$(OBJEXT) dirax.$(OBJEXT) $(am__objects_1) + sfac.$(OBJEXT) dirax.$(OBJEXT) reflections.$(OBJEXT) \ + $(am__objects_1) indexamajig_OBJECTS = $(am_indexamajig_OBJECTS) indexamajig_DEPENDENCIES = am__pattern_sim_SOURCES_DIST = pattern_sim.c diffraction.c utils.c \ @@ -231,7 +232,7 @@ process_hkl_SOURCES = process_hkl.c sfac.c statistics.c cell.c utils.c \ process_hkl_LDADD = @LIBS@ indexamajig_SOURCES = indexamajig.c hdf5-file.c utils.c cell.c image.c \ intensities.c peaks.c index.c filters.c diffraction.c \ - detector.c sfac.c dirax.c $(am__append_3) + detector.c sfac.c dirax.c reflections.c $(am__append_3) indexamajig_LDADD = @LIBS@ @HAVE_GTK_TRUE@hdfsee_SOURCES = hdfsee.c displaywindow.c render.c hdf5-file.c utils.c image.c \ @HAVE_GTK_TRUE@ filters.c @@ -615,3 +615,45 @@ double resolution(UnitCell *cell, signed int h, signed int k, signed int l) return oneoverd / 2; } + + +UnitCell *load_cell_from_pdb(const char *filename) +{ + FILE *fh; + char *rval; + UnitCell *cell = NULL; + + fh = fopen(filename, "r"); + + do { + + char line[1024]; + + rval = fgets(line, 1023, fh); + + if ( strncmp(line, "CRYST1", 6) == 0 ) { + + float a, b, c, al, be, ga; + int r; + + r = sscanf(line+7, "%f %f %f %f %f %f", + &a, &b, &c, &al, &be, &ga); + if ( r != 6 ) { + ERROR("Couldn't understand CRYST1 line\n"); + return NULL; + } + + cell = cell_new_from_parameters(a*1e-10, + b*1e-10, c*1e-10, + deg2rad(al), + deg2rad(be), + deg2rad(ga)); + + } + + } while ( rval != NULL ); + + fclose(fh); + + return cell; +} @@ -77,4 +77,6 @@ extern void cell_print(UnitCell *cell); extern UnitCell *match_cell(UnitCell *cell, UnitCell *template, int verbose); +extern UnitCell *load_cell_from_pdb(const char *filename); + #endif /* CELL_H */ diff --git a/src/diffraction-gpu.c b/src/diffraction-gpu.c index cc7b9e5b..741f7b69 100644 --- a/src/diffraction-gpu.c +++ b/src/diffraction-gpu.c @@ -107,7 +107,7 @@ static void check_sinc_lut(struct gpu_context *gctx, int n) void get_diffraction_gpu(struct gpu_context *gctx, struct image *image, - int na, int nb, int nc, int no_sfac) + int na, int nb, int nc, UnitCell *ucell) { cl_int err; double ax, ay, az; @@ -130,9 +130,7 @@ void get_diffraction_gpu(struct gpu_context *gctx, struct image *image, return; } - cell_get_cartesian(image->molecule->cell, &ax, &ay, &az, - &bx, &by, &bz, - &cx, &cy, &cz); + cell_get_cartesian(ucell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz); cell[0] = ax; cell[1] = ay; cell[2] = az; cell[3] = bx; cell[4] = by; cell[5] = bz; cell[6] = cx; cell[7] = cy; cell[8] = cz; diff --git a/src/diffraction-gpu.h b/src/diffraction-gpu.h index b23dab63..a9e46b66 100644 --- a/src/diffraction-gpu.h +++ b/src/diffraction-gpu.h @@ -24,7 +24,7 @@ struct gpu_context; #if HAVE_OPENCL extern void get_diffraction_gpu(struct gpu_context *gctx, struct image *image, - int na, int nb, int nc); + int na, int nb, int nc, UnitCell *ucell); extern struct gpu_context *setup_gpu(int no_sfac, struct image *image, double *intensities); extern void cleanup_gpu(struct gpu_context *gctx); @@ -32,7 +32,7 @@ extern void cleanup_gpu(struct gpu_context *gctx); #else static void get_diffraction_gpu(struct gpu_context *gctx, struct image *image, - int na, int nb, int nc) + int na, int nb, int nc, UnitCell *ucell) { /* Do nothing */ ERROR("This copy of CrystFEL was not compiled with OpenCL support.\n"); diff --git a/src/diffraction.c b/src/diffraction.c index 01624be8..c47920f5 100644 --- a/src/diffraction.c +++ b/src/diffraction.c @@ -171,7 +171,7 @@ struct rvec get_q(struct image *image, unsigned int xs, unsigned int ys, void get_diffraction(struct image *image, int na, int nb, int nc, - double *intensities, int do_water) + double *intensities, UnitCell *cell, int do_water) { unsigned int xs, ys; double ax, ay, az; @@ -179,11 +179,7 @@ void get_diffraction(struct image *image, int na, int nb, int nc, double cx, cy, cz; float k, klow, bwstep; - if ( image->molecule == NULL ) return; - - cell_get_cartesian(image->molecule->cell, &ax, &ay, &az, - &bx, &by, &bz, - &cx, &cy, &cz); + cell_get_cartesian(cell, &ax, &ay, &az, &bx, &by, &bz, &cx, &cy, &cz); /* Allocate (and zero) the "diffraction array" */ image->data = calloc(image->width * image->height, sizeof(float)); diff --git a/src/diffraction.h b/src/diffraction.h index dd0795fb..58ae10c3 100644 --- a/src/diffraction.h +++ b/src/diffraction.h @@ -20,7 +20,7 @@ #include "cell.h" extern void get_diffraction(struct image *image, int na, int nb, int nc, - double *intensities, int do_water); + double *intensities, UnitCell *cell, int do_water); extern struct rvec get_q(struct image *image, unsigned int xs, unsigned int ys, unsigned int sampling, float *ttp, float k); diff --git a/src/get_hkl.c b/src/get_hkl.c index 29ff0786..5e5331b4 100644 --- a/src/get_hkl.c +++ b/src/get_hkl.c @@ -110,8 +110,6 @@ int main(int argc, char *argv[]) char *output = NULL; unsigned int *counts; signed int h, k, l; - FILE *fh; - char *rval; /* Long options */ const struct option longopts[] = { diff --git a/src/image.h b/src/image.h index 5f2cf3d0..2885b734 100644 --- a/src/image.h +++ b/src/image.h @@ -72,7 +72,6 @@ struct image { float *data; double *twotheta; - struct molecule *molecule; UnitCell *indexed_cell; UnitCell *candidate_cells[MAX_CELL_CANDIDATES]; int ncells; diff --git a/src/index.c b/src/index.c index daa16b27..f2b8ffb9 100644 --- a/src/index.c +++ b/src/index.c @@ -57,8 +57,8 @@ static void write_drx(struct image *image) } -void index_pattern(struct image *image, IndexingMethod indm, int no_match, - int verbose) +void index_pattern(struct image *image, UnitCell *cell, IndexingMethod indm, + int no_match, int verbose) { int i; int nc = 0; @@ -115,8 +115,7 @@ void index_pattern(struct image *image, IndexingMethod indm, int no_match, STATUS("--------------------\n"); } - new_cell = match_cell(image->candidate_cells[i], - image->molecule->cell, verbose); + new_cell = match_cell(image->candidate_cells[i], cell, verbose); image->indexed_cell = new_cell; if ( new_cell != NULL ) { STATUS("Matched on attempt %i.\n", i); diff --git a/src/index.h b/src/index.h index 26ce1c44..b30db837 100644 --- a/src/index.h +++ b/src/index.h @@ -18,6 +18,9 @@ #endif +#include "cell.h" + + typedef enum { INDEXING_NONE, INDEXING_DIRAX, @@ -25,7 +28,7 @@ typedef enum { } IndexingMethod; -extern void index_pattern(struct image *image, IndexingMethod indm, - int no_match, int verbose); +extern void index_pattern(struct image *image, UnitCell *cell, + IndexingMethod indm, int no_match, int verbose); #endif /* INDEX_H */ diff --git a/src/indexamajig.c b/src/indexamajig.c index 5a47254f..e0871fe7 100644 --- a/src/indexamajig.c +++ b/src/indexamajig.c @@ -32,6 +32,7 @@ #include "detector.h" #include "sfac.h" #include "filters.h" +#include "reflections.h" static void show_help(const char *s) @@ -77,6 +78,8 @@ static void show_help(const char *s) " --no-match Don't attempt to match the indexed cell to the\n" " model, just proceed with the one generated by the\n" " auto-indexing procedure.\n" +" --intensities=<file> Specify file containing reflection intensities\n" +" to use.\n" ); } @@ -138,30 +141,22 @@ static struct image *get_simage(struct image *template, int alternate) image->lambda = ph_en_to_lambda(eV_to_J(1.8e3)); - image->molecule = copy_molecule(template->molecule); - free(image->molecule->cell); - image->molecule->cell = cell_new_from_cell(template->indexed_cell); - return image; } -static void simulate_and_write(struct image *simage, - struct gpu_context **gctx) +static void simulate_and_write(struct image *simage, struct gpu_context **gctx, + double *intensities, UnitCell *cell) { /* Set up GPU if necessary */ if ( (gctx != NULL) && (*gctx == NULL) ) { - *gctx = setup_gpu(0, simage, simage->molecule); + *gctx = setup_gpu(0, simage, intensities); } if ( (gctx != NULL) && (*gctx != NULL) ) { - get_diffraction_gpu(*gctx, simage, 24, 24, 40); + get_diffraction_gpu(*gctx, simage, 24, 24, 40, cell); } else { - get_diffraction(simage, 24, 24, 40, 0, 0); - } - if ( simage->molecule == NULL ) { - ERROR("Couldn't open molecule.pdb\n"); - return; + get_diffraction(simage, 24, 24, 40, intensities, cell, 0); } record_image(simage, 0); @@ -193,6 +188,9 @@ int main(int argc, char *argv[]) IndexingMethod indm; char *indm_str = NULL; struct image image; + UnitCell *cell; + double *intensities = NULL; + char *intfile = NULL; /* Long options */ const struct option longopts[] = { @@ -210,6 +208,7 @@ int main(int argc, char *argv[]) {"no-match", 0, &config_nomatch, 1}, {"verbose", 0, &config_verbose, 1}, {"alternate", 0, &config_alternate, 1}, + {"intensities", 0, NULL, 'q'}, {0, 0, NULL, 0} }; @@ -232,6 +231,11 @@ int main(int argc, char *argv[]) break; } + case 'q' : { + intfile = strdup(optarg); + break; + } + case 0 : { break; } @@ -257,6 +261,12 @@ int main(int argc, char *argv[]) return 1; } + if ( intfile != NULL ) { + intensities = read_reflections(intfile); + } else { + intensities = NULL; + } + if ( indm_str == NULL ) { STATUS("You didn't specify an indexing method, so I won't" " try to index anything.\n" @@ -273,10 +283,9 @@ int main(int argc, char *argv[]) } free(indm_str); - image.molecule = load_molecule(); - if ( image.molecule == NULL ) { - ERROR("You must provide molecule.pdb in the working " - "directory.\n"); + cell = load_cell_from_pdb("molecule.pdb"); + if ( cell == NULL ) { + ERROR("Couldn't read unit cell (from molecule.pdb)\n"); return 1; } @@ -354,7 +363,7 @@ int main(int argc, char *argv[]) /* Calculate orientation matrix (by magic) */ if ( config_writedrx || (indm != INDEXING_NONE) ) { - index_pattern(&image, indm, config_nomatch, + index_pattern(&image, cell, indm, config_nomatch, config_verbose); } @@ -376,16 +385,17 @@ int main(int argc, char *argv[]) /* Simulate if requested */ if ( config_simulate ) { if ( config_gpu ) { - simulate_and_write(simage, &gctx); + simulate_and_write(simage, &gctx, intensities, + cell); } else { - simulate_and_write(simage, NULL); + simulate_and_write(simage, NULL, intensities, + cell); } } /* Finished with alternate image */ if ( simage->twotheta != NULL ) free(simage->twotheta); if ( simage->data != NULL ) free(simage->data); - free_molecule(simage->molecule); free(simage); /* Only free cell if found */ @@ -402,7 +412,6 @@ done: } while ( rval != NULL ); fclose(fh); - free_molecule(image.molecule); STATUS("There were %i images.\n", n_images); STATUS("%i hits were found.\n", n_hits); diff --git a/src/pattern_sim.c b/src/pattern_sim.c index 4fd7da66..2ff12226 100644 --- a/src/pattern_sim.c +++ b/src/pattern_sim.c @@ -162,6 +162,7 @@ int main(int argc, char *argv[]) int number = 1; /* Number used for filename of image */ int n_images = 1; /* Generate one image by default */ int done = 0; + UnitCell *cell; /* Long options */ const struct option longopts[] = { @@ -240,7 +241,7 @@ int main(int argc, char *argv[]) image.width = 1024; image.height = 1024; image.lambda = ph_en_to_lambda(eV_to_J(1790.0)); /* Wavelength */ - image.molecule = load_molecule(); + cell = load_cell_from_pdb("molecule.pdb"); #include "geometry-lcls.tmp" @@ -281,7 +282,7 @@ int main(int argc, char *argv[]) image.data = NULL; image.twotheta = NULL; - cell_get_parameters(image.molecule->cell, &a, &b, &c, &d, &d, &d); + cell_get_parameters(cell, &a, &b, &c, &d, &d, &d); STATUS("Particle size = %i x %i x %i (=%5.2f x %5.2f x %5.2f nm)\n", na, nb, nc, na*a/1.0e-9, nb*b/1.0e-9, nc*c/1.0e-9); @@ -290,15 +291,11 @@ int main(int argc, char *argv[]) gctx = setup_gpu(config_nosfac, &image, intensities); } - get_diffraction_gpu(gctx, &image, na, nb, nc); + get_diffraction_gpu(gctx, &image, na, nb, nc, cell); } else { - get_diffraction(&image, na, nb, nc, intensities, + get_diffraction(&image, na, nb, nc, intensities, cell, !config_nowater); } - if ( image.molecule == NULL ) { - ERROR("Couldn't open molecule.pdb\n"); - return 1; - } if ( image.data == NULL ) { ERROR("Diffraction calculation failed.\n"); goto skip; @@ -307,7 +304,7 @@ int main(int argc, char *argv[]) record_image(&image, !config_nonoise); if ( config_nearbragg ) { - output_intensities(&image, image.molecule->cell); + output_intensities(&image, cell); } if ( config_powder ) { @@ -359,8 +356,8 @@ skip: free(image.det.panels); free(powder); - free(image.molecule->reflections); - free(image.molecule); + free(cell); + free(intensities); return 0; } |