Commit bdd2a7b4 authored by Thomas White's avatar Thomas White
Browse files

Allow molecule_factor() to tell when it doesn't know the intensity

parent a6e064b1
......@@ -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;
......
......@@ -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);
......
......@@ -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);
......
......@@ -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);
}
}
......
......@@ -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");
......
......@@ -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;
......
......@@ -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 );
......
......@@ -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);
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment