diff options
| -rw-r--r-- | src/detector.c | 2 | ||||
| -rw-r--r-- | src/utils.c | 22 | ||||
| -rw-r--r-- | src/utils.h | 1 |
3 files changed, 25 insertions, 0 deletions
diff --git a/src/detector.c b/src/detector.c index 1a7d36ee..fae37126 100644 --- a/src/detector.c +++ b/src/detector.c @@ -177,6 +177,8 @@ void record_image(struct image *image) counts = intensity * ph_per_e * sa; //printf("%e counts\n", counts); + counts += poisson_noise(counts); + image->hdr[x + image->width*y] = counts; } diff --git a/src/utils.c b/src/utils.c index d5e1b874..faf78502 100644 --- a/src/utils.c +++ b/src/utils.c @@ -71,3 +71,25 @@ void progress_bar(int val, int total) fflush(stdout); } + + +double poisson_noise(double expected) +{ + double L; + int k; + double p = 1.0; + + L = exp(-expected); + + do { + + double r; + + k++; + r = (double)random()/RAND_MAX; + p *= r; + + } while ( p > L ); + + return k - 1; +} diff --git a/src/utils.h b/src/utils.h index 9e6458cf..54cfdcfb 100644 --- a/src/utils.h +++ b/src/utils.h @@ -55,6 +55,7 @@ extern double angle_between(double x1, double y1, double z1, extern size_t skipspace(const char *s); extern void chomp(char *s); extern void progress_bar(int val, int total); +extern double poisson_noise(double expected); /* Keep these ones inline, to avoid function call overhead */ static inline double distance(double x1, double y1, double x2, double y2) |