aboutsummaryrefslogtreecommitdiff
path: root/src/detector.c
blob: 13d5a20b9dcf79922b35ae5d8cc2f95f438963c3 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
/*
 * detector.c
 *
 * Detector properties
 *
 * (c) 2007-2009 Thomas White <thomas.white@desy.de>
 *
 * pattern_sim - Simulate diffraction patterns from small crystals
 *
 */


#include <stdlib.h>
#include <math.h>
#include <stdio.h>

#include "image.h"
#include "utils.h"


/* Pulse energy density in J/m^2 */
#define PULSE_ENERGY_DENSITY (30.0e7)


/* Detector's quantum efficiency */
#define DQE (0.8)


static uint16_t *bloom(double *hdr, int width, int height)
{
	int x, y;
	uint16_t *data;
	
	data = malloc(width * height * sizeof(uint16_t));
	
	for ( x=0; x<width; x++ ) {
	for ( y=0; y<height; y++ ) {
	
		double hdval;
		
		hdval = hdr[x + width*y] * DQE;
		
		
	
	}
	}
	
	return data;
}


void record_image(struct image *image)
{
	int x, y;
	double ph_per_e;
	double twotheta_max, np, sa_per_pixel;

	/* How many photons are scattered per electron? */
	ph_per_e = PULSE_ENERGY_DENSITY * pow(THOMSON_LENGTH, 2.0)
	          / image->xray_energy;
	printf("%e photons are scattered per electron\n", ph_per_e);

	twotheta_max = image->twotheta[0];
	np = sqrt(pow(image->x_centre, 2.0) + pow(image->y_centre, 2.0));
	sa_per_pixel = pow(2.0 * twotheta_max / np, 2.0);
	printf("sa per pixel=%e\n", sa_per_pixel);

	image->hdr = malloc(image->width * image->height * sizeof(double));

	for ( x=0; x<image->width; x++ ) {
	for ( y=0; y<image->height; y++ ) {

		double counts;
		double intensity;
		double sa;
		double complex val;

		val = image->sfacs[x + image->width*y];
		intensity = val * conj(val);

		/* What solid angle is subtended by this pixel? */
		sa = sa_per_pixel * cos(image->twotheta[x + image->width*y]);

		counts = intensity * ph_per_e * sa;

		image->hdr[x + image->width*y] = counts;

	}
	}
	
	image->data = bloom(image->hdr, image->width, image->height);
}