Commit 3849a969 authored by Thomas White's avatar Thomas White
Browse files

Add new program: get_hkl, for generating ideal intensity lists

parent ba72cc21
......@@ -11,4 +11,5 @@ src/*.o
src/.deps
src/pattern_sim
src/process_hkl
src/get_hkl
*~
EXTRA_DIST = configure src/cell.h src/hdf5-file.h src/image.h src/relrod.h \
src/utils.h src/diffraction.h src/detector.h src/ewald.h \
src/sfac.h src/intensities.h
src/sfac.h src/intensities.h src/reflections.h
SUBDIRS = src
bin_PROGRAMS = pattern_sim process_hkl
bin_PROGRAMS = pattern_sim process_hkl get_hkl
AM_CFLAGS = -Wall
pattern_sim_SOURCES = pattern_sim.c diffraction.c utils.c image.c cell.c \
hdf5-file.c ewald.c detector.c sfac.c intensities.c
hdf5-file.c ewald.c detector.c sfac.c intensities.c \
reflections.c
pattern_sim_LDADD = @LIBS@
process_hkl_SOURCES = process_hkl.c sfac.c statistics.c cell.c utils.c
process_hkl_SOURCES = process_hkl.c sfac.c statistics.c cell.c utils.c \
reflections.c
process_hkl_LDADD = @LIBS@
get_hkl_SOURCES = get_hkl.c sfac.c cell.c utils.c reflections.c
get_hkl_LDADD = @LIBS@
/*
* get_hkl.c
*
* Small program to write out a list of h,k,l,I values given a structure
*
* (c) 2006-2009 Thomas White <thomas.white@desy.de>
*
* Part of CrystFEL - crystallography with a FEL
*
*/
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif
#include <stdarg.h>
#include <stdlib.h>
#include <stdio.h>
#include <string.h>
#include <unistd.h>
#include <getopt.h>
#include "utils.h"
#include "sfac.h"
#include "reflections.h"
static void show_help(const char *s)
{
printf("Syntax: %s [options]\n\n", s);
printf(
"Write idealised intensity lists.\n"
"\n"
" -h, --help Display this help message.\n");
}
int main(int argc, char *argv[])
{
int c;
double *ref;
struct molecule *mol;
/* Long options */
const struct option longopts[] = {
{"help", 0, NULL, 'h'},
{0, 0, NULL, 0}
};
/* Short options */
while ((c = getopt_long(argc, argv, "hi:e:r", longopts, NULL)) != -1) {
switch (c) {
case 'h' : {
show_help(argv[0]);
return 0;
}
case 0 : {
break;
}
default : {
return 1;
}
}
}
mol = load_molecule();
get_reflections_cached(mol, eV_to_J(2.0e3));
ref = ideal_intensities(mol->reflections);
write_reflections("results/ideal-reflections.hkl", NULL, ref, 1,
mol->cell);
return 0;
}
......@@ -24,6 +24,7 @@
#include "utils.h"
#include "statistics.h"
#include "sfac.h"
#include "reflections.h"
/* Number of divisions for R vs |q| graphs */
......@@ -124,72 +125,6 @@ static void write_RvsQ(const char *name, double *ref, double *trueref,
}
static void write_reflections(const char *filename, unsigned int *counts,
double *ref, int zone_axis, UnitCell *cell)
{
FILE *fh;
signed int h, k, l;
fh = fopen(filename, "w");
/* Write spacings and angle if zone axis pattern */
if ( zone_axis ) {
double a, b, c, alpha, beta, gamma;
cell_get_parameters(cell, &a, &b, &c, &alpha, &beta, &gamma);
fprintf(fh, "a %5.3f nm\n", a*1e9);
fprintf(fh, "b %5.3f nm\n", b*1e9);
fprintf(fh, "angle %5.3f deg\n", rad2deg(gamma));
fprintf(fh, "scale 10\n");
}
for ( h=-INDMAX; h<INDMAX; h++ ) {
for ( k=-INDMAX; k<INDMAX; k++ ) {
for ( l=-INDMAX; l<INDMAX; l++ ) {
int N;
double F;
if ( counts ) {
N = lookup_count(counts, h, k, l);
if ( N == 0 ) continue;
} else {
N = 1;
}
F = lookup_intensity(ref, h, k, l) / N;
if ( zone_axis && (l != 0) ) continue;
fprintf(fh, "%3i %3i %3i %f\n", h, k, l, F);
}
}
}
fclose(fh);
}
static double *ideal_intensities(double complex *sfac)
{
double *ref;
signed int h, k, l;
ref = new_list_intensity();
/* Generate ideal reflections from complex structure factors */
for ( h=-INDMAX; h<=INDMAX; h++ ) {
for ( k=-INDMAX; k<=INDMAX; k++ ) {
for ( l=-INDMAX; l<=INDMAX; l++ ) {
double complex F = lookup_sfac(sfac, h, k, l);
double intensity = pow(cabs(F), 2.0);
set_intensity(ref, h, k, l, intensity);
}
}
}
return ref;
}
static void process_reflections(double *ref, double *trueref,
unsigned int *counts, unsigned int n_patterns,
UnitCell *cell, int do_rvsq, int do_zoneaxis)
......@@ -315,9 +250,6 @@ int main(int argc, char *argv[])
counts = new_list_count();
trueref = ideal_intensities(mol->reflections);
write_reflections("results/ideal-reflections.hkl", NULL, trueref, 1,
mol->cell);
if ( strcmp(filename, "-") == 0 ) {
fh = stdin;
} else {
......
/*
* reflections.c
*
* Utilities for handling reflections
*
* (c) 2007-2009 Thomas White <thomas.white@desy.de>
*
* Part of CrystFEL - crystallography with a FEL
*
*/
#include <stdlib.h>
#include <math.h>
#include <stdio.h>
#include <complex.h>
#include <string.h>
#include "utils.h"
#include "cell.h"
#include "reflections.h"
void write_reflections(const char *filename, unsigned int *counts,
double *ref, int zone_axis, UnitCell *cell)
{
FILE *fh;
signed int h, k, l;
fh = fopen(filename, "w");
/* Write spacings and angle if zone axis pattern */
if ( zone_axis ) {
double a, b, c, alpha, beta, gamma;
cell_get_parameters(cell, &a, &b, &c, &alpha, &beta, &gamma);
fprintf(fh, "a %5.3f nm\n", a*1e9);
fprintf(fh, "b %5.3f nm\n", b*1e9);
fprintf(fh, "angle %5.3f deg\n", rad2deg(gamma));
fprintf(fh, "scale 10\n");
}
for ( h=-INDMAX; h<INDMAX; h++ ) {
for ( k=-INDMAX; k<INDMAX; k++ ) {
for ( l=-INDMAX; l<INDMAX; l++ ) {
int N;
double F;
if ( counts ) {
N = lookup_count(counts, h, k, l);
if ( N == 0 ) continue;
} else {
N = 1;
}
F = lookup_intensity(ref, h, k, l) / N;
if ( zone_axis && (l != 0) ) continue;
fprintf(fh, "%3i %3i %3i %f\n", h, k, l, F);
}
}
}
fclose(fh);
}
double *ideal_intensities(double complex *sfac)
{
double *ref;
signed int h, k, l;
ref = new_list_intensity();
/* Generate ideal reflections from complex structure factors */
for ( h=-INDMAX; h<=INDMAX; h++ ) {
for ( k=-INDMAX; k<=INDMAX; k++ ) {
for ( l=-INDMAX; l<=INDMAX; l++ ) {
double complex F = lookup_sfac(sfac, h, k, l);
double intensity = pow(cabs(F), 2.0);
set_intensity(ref, h, k, l, intensity);
}
}
}
return ref;
}
/*
* reflections.h
*
* Utilities for handling reflections
*
* (c) 2007-2009 Thomas White <thomas.white@desy.de>
*
* Part of CrystFEL - crystallography with a FEL
*
*/
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif
#ifndef REFLECTIONS_H
#define REFLECTIONS_H
#include "cell.h"
extern void write_reflections(const char *filename, unsigned int *counts,
double *ref, int zone_axis, UnitCell *cell);
extern double *ideal_intensities(double complex *sfac);
#endif /* REFLECTIONS_H */
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