aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorThomas White <taw@bitwiz.org.uk>2010-03-27 13:12:51 +0100
committerThomas White <taw@bitwiz.org.uk>2010-03-27 13:12:51 +0100
commitbdd2a7b4a8a5456a023d54be6fd5d6a011314c7e (patch)
tree3d13a13a57e724ed94505b3b919d40e18c095cc6
parenta6e064b18ed513d341d8793acbc485d7bc291827 (diff)
Allow molecule_factor() to tell when it doesn't know the intensity
-rw-r--r--src/compare_hkl.c4
-rw-r--r--src/diffraction.c25
-rw-r--r--src/diffraction.h3
-rw-r--r--src/indexamajig.c14
-rw-r--r--src/pattern_sim.c9
-rw-r--r--src/process_hkl.c2
-rw-r--r--src/reflections.c3
-rw-r--r--src/reflections.h2
8 files changed, 40 insertions, 22 deletions
diff --git a/src/compare_hkl.c b/src/compare_hkl.c
index 68920a55..50ca6282 100644
--- a/src/compare_hkl.c
+++ b/src/compare_hkl.c
@@ -98,12 +98,12 @@ int main(int argc, char *argv[])
}
cell = load_cell_from_pdb("molecule.pdb");
- ref1 = read_reflections(afile);
+ ref1 = read_reflections(afile, NULL);
if ( ref1 == NULL ) {
ERROR("Couldn't open file '%s'\n", afile);
return 1;
}
- ref2 = read_reflections(bfile);
+ ref2 = read_reflections(bfile, NULL);
if ( ref2 == NULL ) {
ERROR("Couldn't open file '%s'\n", bfile);
return 1;
diff --git a/src/diffraction.c b/src/diffraction.c
index 5ecf04f5..675aa148 100644
--- a/src/diffraction.c
+++ b/src/diffraction.c
@@ -70,10 +70,11 @@ static double lattice_factor(struct rvec q, double ax, double ay, double az,
/* Look up the structure factor for the nearest Bragg condition */
-static double molecule_factor(double *intensities, struct rvec q,
- double ax, double ay, double az,
- double bx, double by, double bz,
- double cx, double cy, double cz)
+static double molecule_factor(double *intensities, unsigned int *counts,
+ struct rvec q,
+ double ax, double ay, double az,
+ double bx, double by, double bz,
+ double cx, double cy, double cz)
{
double hd, kd, ld;
signed int h, k, l;
@@ -86,8 +87,13 @@ static double molecule_factor(double *intensities, struct rvec q,
k = (signed int)rint(kd);
l = (signed int)rint(ld);
- r = lookup_intensity(intensities, h, k, l);
+ if ( lookup_count(counts, h, k, l) == 0 ) {
+ ERROR("Needed intensity for %i %i %i, but don't have it.\n",
+ h, k, l);
+ return 1.0e20;
+ }
+ r = lookup_intensity(intensities, h, k, l);
return r;
}
@@ -171,7 +177,8 @@ 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, UnitCell *cell, int do_water)
+ double *intensities, unsigned int *counts, UnitCell *cell,
+ int do_water)
{
unsigned int xs, ys;
double ax, ay, az;
@@ -225,8 +232,10 @@ void get_diffraction(struct image *image, int na, int nb, int nc,
if ( intensities == NULL ) {
I_molecule = 1.0e10;
} else {
- I_molecule = molecule_factor(intensities, q,
- ax,ay,az,bx,by,bz,cx,cy,cz);
+ I_molecule = molecule_factor(intensities,
+ counts, q,
+ ax,ay,az,
+ bx,by,bz,cx,cy,cz);
}
I_lattice = pow(f_lattice, 2.0);
diff --git a/src/diffraction.h b/src/diffraction.h
index 58ae10c3..c317060f 100644
--- a/src/diffraction.h
+++ b/src/diffraction.h
@@ -20,7 +20,8 @@
#include "cell.h"
extern void get_diffraction(struct image *image, int na, int nb, int nc,
- double *intensities, UnitCell *cell, int do_water);
+ double *intensities, unsigned int *counts,
+ 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/indexamajig.c b/src/indexamajig.c
index c99a2b8b..c24b0e6b 100644
--- a/src/indexamajig.c
+++ b/src/indexamajig.c
@@ -146,7 +146,8 @@ static struct image *get_simage(struct image *template, int alternate)
static void simulate_and_write(struct image *simage, struct gpu_context **gctx,
- double *intensities, UnitCell *cell)
+ double *intensities, unsigned int *counts,
+ UnitCell *cell)
{
/* Set up GPU if necessary */
if ( (gctx != NULL) && (*gctx == NULL) ) {
@@ -156,7 +157,8 @@ static void simulate_and_write(struct image *simage, struct gpu_context **gctx,
if ( (gctx != NULL) && (*gctx != NULL) ) {
get_diffraction_gpu(*gctx, simage, 24, 24, 40, cell);
} else {
- get_diffraction(simage, 24, 24, 40, intensities, cell, 0);
+ get_diffraction(simage, 24, 24, 40,
+ intensities, counts, cell, 0);
}
record_image(simage, 0);
@@ -191,6 +193,7 @@ int main(int argc, char *argv[])
UnitCell *cell;
double *intensities = NULL;
char *intfile = NULL;
+ unsigned int *counts = NULL;
/* Long options */
const struct option longopts[] = {
@@ -262,9 +265,10 @@ int main(int argc, char *argv[])
}
if ( intfile != NULL ) {
- intensities = read_reflections(intfile);
+ intensities = read_reflections(intfile, counts);
} else {
intensities = NULL;
+ counts = NULL;
}
if ( indm_str == NULL ) {
@@ -386,10 +390,10 @@ int main(int argc, char *argv[])
if ( config_simulate ) {
if ( config_gpu ) {
simulate_and_write(simage, &gctx, intensities,
- cell);
+ counts, cell);
} else {
simulate_and_write(simage, NULL, intensities,
- cell);
+ counts, cell);
}
}
diff --git a/src/pattern_sim.c b/src/pattern_sim.c
index 55b92c0e..3fa1155f 100644
--- a/src/pattern_sim.c
+++ b/src/pattern_sim.c
@@ -163,6 +163,7 @@ int main(int argc, char *argv[])
int n_images = 1; /* Generate one image by default */
int done = 0;
UnitCell *cell;
+ unsigned int *counts;
/* Long options */
const struct option longopts[] = {
@@ -232,8 +233,10 @@ int main(int argc, char *argv[])
STATUS("reflection intensities (with --intensities).\n");
STATUS("I'll simulate a flat intensity distribution.\n");
intensities = NULL;
+ counts = NULL;
} else {
- intensities = read_reflections(intfile);
+ counts = new_list_count();
+ intensities = read_reflections(intfile, counts);
free(intfile);
}
@@ -293,8 +296,8 @@ int main(int argc, char *argv[])
}
get_diffraction_gpu(gctx, &image, na, nb, nc, cell);
} else {
- get_diffraction(&image, na, nb, nc, intensities, cell,
- !config_nowater);
+ get_diffraction(&image, na, nb, nc, intensities, counts,
+ cell, !config_nowater);
}
if ( image.data == NULL ) {
ERROR("Diffraction calculation failed.\n");
diff --git a/src/process_hkl.c b/src/process_hkl.c
index 7b1035bc..8bc76e5e 100644
--- a/src/process_hkl.c
+++ b/src/process_hkl.c
@@ -268,7 +268,7 @@ int main(int argc, char *argv[])
if ( intfile != NULL ) {
STATUS("Comparing against '%s'\n", intfile);
- trueref = read_reflections(intfile);
+ trueref = read_reflections(intfile, NULL);
free(intfile);
} else {
trueref = NULL;
diff --git a/src/reflections.c b/src/reflections.c
index deba0f8d..c943d014 100644
--- a/src/reflections.c
+++ b/src/reflections.c
@@ -84,7 +84,7 @@ void write_reflections(const char *filename, unsigned int *counts,
}
-double *read_reflections(const char *filename)
+double *read_reflections(const char *filename, unsigned int *counts)
{
double *ref;
FILE *fh;
@@ -110,6 +110,7 @@ double *read_reflections(const char *filename)
if ( r != 4 ) continue;
set_intensity(ref, h, k, l, intensity);
+ if ( counts != NULL ) set_count(counts, h, k, l, 1);
} while ( rval != NULL );
diff --git a/src/reflections.h b/src/reflections.h
index 14209c34..8c513a64 100644
--- a/src/reflections.h
+++ b/src/reflections.h
@@ -23,7 +23,7 @@
extern void write_reflections(const char *filename, unsigned int *counts,
double *ref, int zone_axis, UnitCell *cell);
-extern double *read_reflections(const char *filename);
+extern double *read_reflections(const char *filename, unsigned int *counts);
extern double *ideal_intensities(double complex *sfac);