diff options
-rw-r--r-- | src/find_peaks.c | 107 | ||||
-rw-r--r-- | src/find_peaks.h | 4 | ||||
-rw-r--r-- | src/fit.c | 6 | ||||
-rw-r--r-- | src/misc.c | 56 | ||||
-rw-r--r-- | src/misc.h | 3 | ||||
-rw-r--r-- | src/test-find-peaks.c | 58 | ||||
-rw-r--r-- | src/test.c | 86 | ||||
-rw-r--r-- | src/viridis.pal | 287 |
8 files changed, 568 insertions, 39 deletions
diff --git a/src/find_peaks.c b/src/find_peaks.c index 05d8c15..e9030a0 100644 --- a/src/find_peaks.c +++ b/src/find_peaks.c @@ -6,6 +6,7 @@ #include <stdlib.h> /* for malloc() */ #include "optics.h" #include <math.h> /* for exp() */ +#include "misc.h" typedef struct peak { size_t i; @@ -115,14 +116,37 @@ end: return; } -void get_hough_transform(event *ev, double *pos, double *x, double *y, size_t n, size_t m, double *result) +void get_hough_transform(event *ev, double *pos, double *x, double *y, size_t n, size_t m, double *result, double *last, size_t len) { /* Computes the "Hough transform" of the event `ev` and stores it in `result`. */ size_t i, j, k; - double dir[3], pmt_dir[3], cos_theta; + double dir[3], pmt_dir[3], cos_theta2; + int skip; + double *sin_theta, *cos_theta, *sin_phi, *cos_phi; + double wavelength0, n_d2o; - double wavelength0 = 400.0; - double n_d2o = get_index_snoman_d2o(wavelength0); + sin_theta = malloc(sizeof(double)*n); + cos_theta = malloc(sizeof(double)*n); + sin_phi = malloc(sizeof(double)*n); + cos_phi = malloc(sizeof(double)*n); + + wavelength0 = 400.0; + n_d2o = get_index_snoman_d2o(wavelength0); + + /* Precompute sin(theta), cos(theta), sin(phi), and cos(phi) to speed + * things up. */ + for (i = 0; i < n; i++) { + sin_theta[i] = sin(x[i]); + cos_theta[i] = cos(x[i]); + } + + for (i = 0; i < m; i++) { + sin_phi[i] = sin(y[i]); + cos_phi[i] = cos(y[i]); + } + + /* Zero out the result array. */ + for (i = 0; i < n*m; i++) result[i] = 0.0; for (i = 0; i < MAX_PMTS; i++) { if (pmts[i].pmt_type != PMT_NORMAL || !ev->pmt_hits[i].hit) continue; @@ -131,55 +155,82 @@ void get_hough_transform(event *ev, double *pos, double *x, double *y, size_t n, normalize(pmt_dir); + /* Ignore PMT hits within the Cerenkov cone of previously found peaks. */ + skip = 0; + for (j = 0; j < len; j++) { + if (DOT(pmt_dir,last+3*j) > 1/n_d2o) { + skip = 1; + break; + } + } + + if (skip) continue; + for (j = 0; j < n; j++) { for (k = 0; k < m; k++) { - double r = 1+x[j]*x[j]+y[k]*y[k]; - dir[0] = 2*x[j]/r; - dir[1] = 2*y[k]/r; - dir[2] = (-1+x[j]*x[j]+y[k]*y[k])/r; + dir[0] = sin_theta[j]*cos_phi[k]; + dir[1] = sin_theta[j]*sin_phi[k]; + dir[2] = cos_theta[j]; - cos_theta = DOT(pmt_dir,dir); + cos_theta2 = DOT(pmt_dir,dir); - result[j*m + k] += ev->pmt_hits[i].qhs*exp(-fabs(cos_theta-1.0/n_d2o)/0.1); + result[j*n + k] += ev->pmt_hits[i].qhs*exp(-fabs(cos_theta2-1.0/n_d2o)/0.1); } } } + + free(sin_theta); + free(cos_theta); + free(sin_phi); + free(cos_phi); } -void find_peaks(event *ev, double *pos, size_t n, size_t m, double *peak_theta, double *peak_phi, size_t *npeaks, size_t max_peaks, double threshold) +void find_peaks(event *ev, double *pos, size_t n, size_t m, double *peak_theta, double *peak_phi, size_t *npeaks, size_t max_peaks) { + /* Finds rings in the event `ev` by looking for peaks in the Hough + * transform. The directions of potential particle rings are stored in the + * arrays `peak_theta` and `peak_phi`. The number of rings found are stored + * in the variable `npeaks`. + * + * This function first computes the Hough transform of the event `ev` and + * looks for the highest peak. The next peak is then found by computing the + * Hough transform ignoring any PMTs within the Cerenkov cone of the first + * peak. This process is repeated until `max_peaks` peaks are found. */ size_t i; + double *x, *y, *result, *dir; + size_t *imax, *jmax; + size_t max; - double *x = calloc(n,sizeof(double)); - double *y = calloc(m,sizeof(double)); - double *result = calloc(n*m,sizeof(double)); + x = calloc(n,sizeof(double)); + y = calloc(m,sizeof(double)); + result = malloc(n*m*sizeof(double)); + dir = calloc(max_peaks*3,sizeof(double)); for (i = 0; i < n; i++) { - x[i] = -10 + 20.0*i/(n-1); + x[i] = i*M_PI/(n-1); } for (i = 0; i < m; i++) { - y[i] = -10 + 20.0*i/(m-1); + y[i] = i*2*M_PI/(m-1); } - get_hough_transform(ev,pos,x,y,n,m,result); + imax = calloc(max_peaks,sizeof(size_t)); + jmax = calloc(max_peaks,sizeof(size_t)); - size_t *imax = calloc(max_peaks,sizeof(size_t)); - size_t *jmax = calloc(max_peaks,sizeof(size_t)); - - find_peaks_array(result,n,m,imax,jmax,npeaks,max_peaks,threshold); - - for (i = 0; i < *npeaks; i++) { - double R, theta; - R = sqrt(x[imax[i]]*x[imax[i]]+y[jmax[i]]*y[jmax[i]]); - theta = atan2(y[jmax[i]],x[imax[i]]); - peak_theta[i] = 2*atan(1/R); - peak_phi[i] = theta; + for (i = 0; i < max_peaks; i++) { + get_hough_transform(ev,pos,x,y,n,m,result,dir,i); + max = argmax(result,n*m); + peak_theta[i] = x[max/m]; + peak_phi[i] = y[max % m]; + get_dir(dir+i*3,peak_theta[i],peak_phi[i]); } + *npeaks = max_peaks; + free(imax); free(jmax); free(x); free(y); free(result); + free(dir); } diff --git a/src/find_peaks.h b/src/find_peaks.h index 2abf6fe..7ca798b 100644 --- a/src/find_peaks.h +++ b/src/find_peaks.h @@ -5,7 +5,7 @@ #include <stddef.h> /* for size_t */ void find_peaks_array(double *x, size_t n, size_t m, size_t *imax, size_t *jmax, size_t *npeaks, size_t max_peaks, double threshold); -void get_hough_transform(event *ev, double *pos, double *x, double *y, size_t n, size_t m, double *result); -void find_peaks(event *ev, double *pos, size_t n, size_t m, double *peak_theta, double *peak_phi, size_t *npeaks, size_t max_peaks, double threshold); +void get_hough_transform(event *ev, double *pos, double *x, double *y, size_t n, size_t m, double *result, double *last, size_t len); +void find_peaks(event *ev, double *pos, size_t n, size_t m, double *peak_theta, double *peak_phi, size_t *npeaks, size_t max_peaks); #endif @@ -5296,7 +5296,7 @@ int fit_event2(event *ev, double *xopt, double *fmin, int *id, size_t n, double ub[3] = GTVALID; /* Find the peaks in the Hough transform of the event. */ - find_peaks(ev, pos, 100, 100, peak_theta, peak_phi, &npeaks, LEN(peak_theta), 0.5); + find_peaks(ev, pos, 100, 100, peak_theta, peak_phi, &npeaks, LEN(peak_theta)); result = malloc(sizeof(size_t)*n*ipow(npeaks,n)); @@ -5321,7 +5321,7 @@ int fit_event2(event *ev, double *xopt, double *fmin, int *id, size_t n, double fpars.hit_only = 0; fpars.print = 0; fpars.charge_only = 0; - nlopt_set_ftol_abs(opt, 1e-2); + nlopt_set_ftol_abs(opt, 1e-5); nlopt_set_maxeval(opt, 10000); fpars.n = n; @@ -5365,7 +5365,7 @@ int fit_event2(event *ev, double *xopt, double *fmin, int *id, size_t n, double if (T0 < Tmin2) T0 = Tmin2; x[4+3*j] = T0; - ss[4+3*j] = x[4+3*j]*0.02; + ss[4+3*j] = x[4+3*j]*0.1; ss[5+3*j] = 0.01; ss[6+3*j] = 0.01; lb[4+3*j] = Tmin; @@ -742,3 +742,59 @@ void combinations_with_replacement(size_t n, size_t r, size_t *result, size_t *l free(tmp); } + +size_t argmax(double *a, size_t n) +{ + /* Returns the index of the maximum of the array `a`. */ + size_t i, index; + double max; + + if (n < 2) return 0; + + index = 0; + max = a[0]; + for (i = 1; i < n; i++) { + if (a[i] > max) { + max = a[i]; + index = i; + } + } + + return index; +} + +size_t argmin(double *a, size_t n) +{ + /* Returns the index of the minimum of the array `a`. */ + size_t i, index; + double min; + + if (n < 2) return 0; + + index = 0; + min = a[0]; + for (i = 1; i < n; i++) { + if (a[i] < min) { + min = a[i]; + index = i; + } + } + + return index; +} + +void get_dir(double *dir, double theta, double phi) +{ + /* Compute the 3-dimensional unit vector `dir` from the polar angle `theta` + * and azimuthal angle `phi`. */ + double sin_theta, cos_theta, sin_phi, cos_phi; + + cos_theta = cos(theta); + sin_theta = sin(theta); + cos_phi = cos(phi); + sin_phi = sin(phi); + + dir[0] = sin_theta*cos_phi; + dir[1] = sin_theta*sin_phi; + dir[2] = cos_theta; +} @@ -31,5 +31,8 @@ size_t ipow(size_t base, size_t exp); void product(size_t n, size_t r, size_t *result); void unique_vertices(int *id, size_t n, size_t npeaks, size_t *result, size_t *nvertices); void combinations_with_replacement(size_t n, size_t r, size_t *result, size_t *len); +size_t argmax(double *a, size_t n); +size_t argmin(double *a, size_t n); +void get_dir(double *dir, double theta, double phi); #endif diff --git a/src/test-find-peaks.c b/src/test-find-peaks.c index c7c2c1d..c82f9ec 100644 --- a/src/test-find-peaks.c +++ b/src/test-find-peaks.c @@ -9,12 +9,13 @@ #include <math.h> /* for M_PI */ #include <errno.h> /* for errno */ #include <string.h> /* for strerror() */ +#include "vector.h" #define EV_RECORD 0x45562020 #define MCTK_RECORD 0x4d43544b #define MCVX_RECORD 0x4d435658 -void plot_3d(double *x, double *y, double *z, size_t n, size_t m) +void plot_hough_transform(double *x, double *y, double *z, size_t n, size_t m) { size_t i, j; @@ -49,6 +50,51 @@ void plot_3d(double *x, double *y, double *z, size_t n, size_t m) } } +void plot_find_peaks(event *ev, double *peak_theta, double *peak_phi, size_t n) +{ + size_t i; + + FILE *pipe = popen("gnuplot -p", "w"); + + if (!pipe) { + fprintf(stderr, "error running gnuplot command: %s\n", strerror(errno)); + exit(1); + } + + fprintf(pipe, "set macros\n"); + fprintf(pipe, "load 'viridis.pal'\n"); + fprintf(pipe, "set title 'Hough Transform'\n"); + /* Not entirely sure what these do, but following the instructions from + * http://lowrank.net/gnuplot/plot3d2-e.html. */ + fprintf(pipe, "set xrange [0:6.28]\n"); + fprintf(pipe, "set yrange [3.14:0]\n"); + fprintf(pipe, "plot '-' u 1:2:3:4 with circles palette fillstyle solid, '-' u 1:2 with circles lc rgb \"red\" lw 2\n"); + + for (i = 0; i < MAX_PMTS; i++) { + if (ev->pmt_hits[i].hit) { + double r = NORM(pmts[i].pos); + double theta = acos(pmts[i].pos[2]/r); + double phi = atan2(pmts[i].pos[1],pmts[i].pos[0]); + phi = (phi < 0) ? phi + 2*M_PI: phi; + fprintf(pipe, "%.10g %.10g %.10g %.10g\n", phi, theta, 0.01, ev->pmt_hits[i].qhs); + } + } + fprintf(pipe,"e\n"); + + for (i = 0; i < n; i++) { + fprintf(pipe, "%.10g %.10g\n", peak_phi[i], peak_theta[i]); + } + fprintf(pipe,"e\n"); + + /* Pause so that you can rotate the 3D graph. */ + fprintf(pipe,"pause mouse keypress\n"); + + if (pclose(pipe)) { + fprintf(stderr, "error closing gnuplot command: %s\n", strerror(errno)); + exit(1); + } +} + int get_event(zebraFile *f, event *ev, bank *b) { /* Read all the PMT banks from the zebra file and update `ev`. @@ -193,16 +239,16 @@ int main(int argc, char **argv) double *result = calloc(n*m,sizeof(double)); for (i = 0; i < n; i++) { - x[i] = -10 + 20.0*i/(n-1); + x[i] = i*M_PI/(n-1); } for (i = 0; i < m; i++) { - y[i] = -10 + 20.0*i/(m-1); + y[i] = i*2*M_PI/(m-1); } - get_hough_transform(&ev,pos,x,y,n,m,result); + get_hough_transform(&ev,pos,x,y,n,m,result,0,0); - find_peaks(&ev,pos,n,m,peak_theta,peak_phi,&npeaks,LEN(peak_theta),0.5); + find_peaks(&ev,pos,n,m,peak_theta,peak_phi,&npeaks,LEN(peak_theta)); printf("gtid %i\n", ev.gtid); for (i = 0; i < npeaks; i++) { @@ -210,7 +256,7 @@ int main(int argc, char **argv) } if (plot) - plot_3d(x,y,result,n,m); + plot_find_peaks(&ev,peak_theta,peak_phi,npeaks); /* Skip reading in the next bank since get_event() already read in * the next bank. */ @@ -25,6 +25,7 @@ #include "sno.h" #include "quad.h" #include "find_peaks.h" +#include <gsl/gsl_sort.h> typedef int testFunction(char *err); @@ -1893,6 +1894,88 @@ int test_combinations_with_replacement(char *err) return 0; } +int test_get_dir(char *err) +{ + /* Tests the get_dir() function. */ + size_t i, j; + double dir[3], expected_dir[3]; + double theta, phi; + + init_genrand(0); + + for (i = 0; i < 100; i++) { + theta = genrand_real2()*M_PI; + phi = genrand_real2()*2*M_PI; + + get_dir(dir,theta,phi); + + expected_dir[0] = sin(theta)*cos(phi); + expected_dir[1] = sin(theta)*sin(phi); + expected_dir[2] = cos(theta); + + for (j = 0; j < 3; j++) { + if (!isclose(dir[j],expected_dir[j],0,1e-9)) { + sprintf(err, "dir[%zu] returned %.2f but expected %.2f", j, dir[j], expected_dir[j]); + return 1; + } + } + } + + return 0; +} + +int test_argmax(char *err) +{ + /* Tests the argmax() function. */ + size_t i, j; + double a[100]; + size_t max, p; + + init_genrand(0); + + for (i = 0; i < 100; i++) { + for (j = 0; j < 100; j++) + a[j] = genrand_real2(); + + max = argmax(a,LEN(a)); + + gsl_sort_largest_index(&p,1,a,1,LEN(a)); + + if (max != p) { + sprintf(err, "argmax() returned %zu but expected %zu", max, p); + return 1; + } + } + + return 0; +} + +int test_argmin(char *err) +{ + /* Tests the argmin() function. */ + size_t i, j; + double a[100]; + size_t min, p; + + init_genrand(0); + + for (i = 0; i < 100; i++) { + for (j = 0; j < 100; j++) + a[j] = genrand_real2(); + + min = argmin(a,LEN(a)); + + gsl_sort_smallest_index(&p,1,a,1,LEN(a)); + + if (min != p) { + sprintf(err, "argmin() returned %zu but expected %zu", min, p); + return 1; + } + } + + return 0; +} + struct tests { testFunction *test; char *name; @@ -1940,6 +2023,9 @@ struct tests { {test_find_peaks_highest, "test_find_peaks_highest"}, {test_find_peaks_sorted, "test_find_peaks_sorted"}, {test_combinations_with_replacement, "test_combinations_with_replacement"}, + {test_get_dir, "test_get_dir"}, + {test_argmax, "test_argmax"}, + {test_argmin, "test_argmin"}, }; int main(int argc, char **argv) diff --git a/src/viridis.pal b/src/viridis.pal new file mode 100644 index 0000000..ca3efc2 --- /dev/null +++ b/src/viridis.pal @@ -0,0 +1,287 @@ +# New matplotlib colormaps by Nathaniel J. Smith, Stefan van der Walt, +# and (in the case of viridis) Eric Firing. +# +# This file and the colormaps in it are released under the CC0 license / +# public domain dedication. We would appreciate credit if you use or +# redistribute these colormaps, but do not impose any legal restrictions. +# +# To the extent possible under law, the persons who associated CC0 with +# mpl-colormaps have waived all copyright and related or neighboring rights +# to mpl-colormaps. +# +# You should have received a copy of the CC0 legalcode along with this +# work. If not, see <http://creativecommons.org/publicdomain/zero/1.0/>. + +#https://github.com/BIDS/colormap/blob/master/colormaps.py + + +# line styles +set style line 1 lt 1 lc rgb '#440154' # dark purple +set style line 2 lt 1 lc rgb '#472c7a' # purple +set style line 3 lt 1 lc rgb '#3b518b' # blue +set style line 4 lt 1 lc rgb '#2c718e' # blue +set style line 5 lt 1 lc rgb '#21908d' # blue-green +set style line 6 lt 1 lc rgb '#27ad81' # green +set style line 7 lt 1 lc rgb '#5cc863' # green +set style line 8 lt 1 lc rgb '#aadc32' # lime green +set style line 9 lt 1 lc rgb '#fde725' # yellow + + +# palette +set palette defined (\ +0 0.267004 0.004874 0.329415,\ +1 0.268510 0.009605 0.335427,\ +2 0.269944 0.014625 0.341379,\ +3 0.271305 0.019942 0.347269,\ +4 0.272594 0.025563 0.353093,\ +5 0.273809 0.031497 0.358853,\ +6 0.274952 0.037752 0.364543,\ +7 0.276022 0.044167 0.370164,\ +8 0.277018 0.050344 0.375715,\ +9 0.277941 0.056324 0.381191,\ +10 0.278791 0.062145 0.386592,\ +11 0.279566 0.067836 0.391917,\ +12 0.280267 0.073417 0.397163,\ +13 0.280894 0.078907 0.402329,\ +14 0.281446 0.084320 0.407414,\ +15 0.281924 0.089666 0.412415,\ +16 0.282327 0.094955 0.417331,\ +17 0.282656 0.100196 0.422160,\ +18 0.282910 0.105393 0.426902,\ +19 0.283091 0.110553 0.431554,\ +20 0.283197 0.115680 0.436115,\ +21 0.283229 0.120777 0.440584,\ +22 0.283187 0.125848 0.444960,\ +23 0.283072 0.130895 0.449241,\ +24 0.282884 0.135920 0.453427,\ +25 0.282623 0.140926 0.457517,\ +26 0.282290 0.145912 0.461510,\ +27 0.281887 0.150881 0.465405,\ +28 0.281412 0.155834 0.469201,\ +29 0.280868 0.160771 0.472899,\ +30 0.280255 0.165693 0.476498,\ +31 0.279574 0.170599 0.479997,\ +32 0.278826 0.175490 0.483397,\ +33 0.278012 0.180367 0.486697,\ +34 0.277134 0.185228 0.489898,\ +35 0.276194 0.190074 0.493001,\ +36 0.275191 0.194905 0.496005,\ +37 0.274128 0.199721 0.498911,\ +38 0.273006 0.204520 0.501721,\ +39 0.271828 0.209303 0.504434,\ +40 0.270595 0.214069 0.507052,\ +41 0.269308 0.218818 0.509577,\ +42 0.267968 0.223549 0.512008,\ +43 0.266580 0.228262 0.514349,\ +44 0.265145 0.232956 0.516599,\ +45 0.263663 0.237631 0.518762,\ +46 0.262138 0.242286 0.520837,\ +47 0.260571 0.246922 0.522828,\ +48 0.258965 0.251537 0.524736,\ +49 0.257322 0.256130 0.526563,\ +50 0.255645 0.260703 0.528312,\ +51 0.253935 0.265254 0.529983,\ +52 0.252194 0.269783 0.531579,\ +53 0.250425 0.274290 0.533103,\ +54 0.248629 0.278775 0.534556,\ +55 0.246811 0.283237 0.535941,\ +56 0.244972 0.287675 0.537260,\ +57 0.243113 0.292092 0.538516,\ +58 0.241237 0.296485 0.539709,\ +59 0.239346 0.300855 0.540844,\ +60 0.237441 0.305202 0.541921,\ +61 0.235526 0.309527 0.542944,\ +62 0.233603 0.313828 0.543914,\ +63 0.231674 0.318106 0.544834,\ +64 0.229739 0.322361 0.545706,\ +65 0.227802 0.326594 0.546532,\ +66 0.225863 0.330805 0.547314,\ +67 0.223925 0.334994 0.548053,\ +68 0.221989 0.339161 0.548752,\ +69 0.220057 0.343307 0.549413,\ +70 0.218130 0.347432 0.550038,\ +71 0.216210 0.351535 0.550627,\ +72 0.214298 0.355619 0.551184,\ +73 0.212395 0.359683 0.551710,\ +74 0.210503 0.363727 0.552206,\ +75 0.208623 0.367752 0.552675,\ +76 0.206756 0.371758 0.553117,\ +77 0.204903 0.375746 0.553533,\ +78 0.203063 0.379716 0.553925,\ +79 0.201239 0.383670 0.554294,\ +80 0.199430 0.387607 0.554642,\ +81 0.197636 0.391528 0.554969,\ +82 0.195860 0.395433 0.555276,\ +83 0.194100 0.399323 0.555565,\ +84 0.192357 0.403199 0.555836,\ +85 0.190631 0.407061 0.556089,\ +86 0.188923 0.410910 0.556326,\ +87 0.187231 0.414746 0.556547,\ +88 0.185556 0.418570 0.556753,\ +89 0.183898 0.422383 0.556944,\ +90 0.182256 0.426184 0.557120,\ +91 0.180629 0.429975 0.557282,\ +92 0.179019 0.433756 0.557430,\ +93 0.177423 0.437527 0.557565,\ +94 0.175841 0.441290 0.557685,\ +95 0.174274 0.445044 0.557792,\ +96 0.172719 0.448791 0.557885,\ +97 0.171176 0.452530 0.557965,\ +98 0.169646 0.456262 0.558030,\ +99 0.168126 0.459988 0.558082,\ +100 0.166617 0.463708 0.558119,\ +101 0.165117 0.467423 0.558141,\ +102 0.163625 0.471133 0.558148,\ +103 0.162142 0.474838 0.558140,\ +104 0.160665 0.478540 0.558115,\ +105 0.159194 0.482237 0.558073,\ +106 0.157729 0.485932 0.558013,\ +107 0.156270 0.489624 0.557936,\ +108 0.154815 0.493313 0.557840,\ +109 0.153364 0.497000 0.557724,\ +110 0.151918 0.500685 0.557587,\ +111 0.150476 0.504369 0.557430,\ +112 0.149039 0.508051 0.557250,\ +113 0.147607 0.511733 0.557049,\ +114 0.146180 0.515413 0.556823,\ +115 0.144759 0.519093 0.556572,\ +116 0.143343 0.522773 0.556295,\ +117 0.141935 0.526453 0.555991,\ +118 0.140536 0.530132 0.555659,\ +119 0.139147 0.533812 0.555298,\ +120 0.137770 0.537492 0.554906,\ +121 0.136408 0.541173 0.554483,\ +122 0.135066 0.544853 0.554029,\ +123 0.133743 0.548535 0.553541,\ +124 0.132444 0.552216 0.553018,\ +125 0.131172 0.555899 0.552459,\ +126 0.129933 0.559582 0.551864,\ +127 0.128729 0.563265 0.551229,\ +128 0.127568 0.566949 0.550556,\ +129 0.126453 0.570633 0.549841,\ +130 0.125394 0.574318 0.549086,\ +131 0.124395 0.578002 0.548287,\ +132 0.123463 0.581687 0.547445,\ +133 0.122606 0.585371 0.546557,\ +134 0.121831 0.589055 0.545623,\ +135 0.121148 0.592739 0.544641,\ +136 0.120565 0.596422 0.543611,\ +137 0.120092 0.600104 0.542530,\ +138 0.119738 0.603785 0.541400,\ +139 0.119512 0.607464 0.540218,\ +140 0.119423 0.611141 0.538982,\ +141 0.119483 0.614817 0.537692,\ +142 0.119699 0.618490 0.536347,\ +143 0.120081 0.622161 0.534946,\ +144 0.120638 0.625828 0.533488,\ +145 0.121380 0.629492 0.531973,\ +146 0.122312 0.633153 0.530398,\ +147 0.123444 0.636809 0.528763,\ +148 0.124780 0.640461 0.527068,\ +149 0.126326 0.644107 0.525311,\ +150 0.128087 0.647749 0.523491,\ +151 0.130067 0.651384 0.521608,\ +152 0.132268 0.655014 0.519661,\ +153 0.134692 0.658636 0.517649,\ +154 0.137339 0.662252 0.515571,\ +155 0.140210 0.665859 0.513427,\ +156 0.143303 0.669459 0.511215,\ +157 0.146616 0.673050 0.508936,\ +158 0.150148 0.676631 0.506589,\ +159 0.153894 0.680203 0.504172,\ +160 0.157851 0.683765 0.501686,\ +161 0.162016 0.687316 0.499129,\ +162 0.166383 0.690856 0.496502,\ +163 0.170948 0.694384 0.493803,\ +164 0.175707 0.697900 0.491033,\ +165 0.180653 0.701402 0.488189,\ +166 0.185783 0.704891 0.485273,\ +167 0.191090 0.708366 0.482284,\ +168 0.196571 0.711827 0.479221,\ +169 0.202219 0.715272 0.476084,\ +170 0.208030 0.718701 0.472873,\ +171 0.214000 0.722114 0.469588,\ +172 0.220124 0.725509 0.466226,\ +173 0.226397 0.728888 0.462789,\ +174 0.232815 0.732247 0.459277,\ +175 0.239374 0.735588 0.455688,\ +176 0.246070 0.738910 0.452024,\ +177 0.252899 0.742211 0.448284,\ +178 0.259857 0.745492 0.444467,\ +179 0.266941 0.748751 0.440573,\ +180 0.274149 0.751988 0.436601,\ +181 0.281477 0.755203 0.432552,\ +182 0.288921 0.758394 0.428426,\ +183 0.296479 0.761561 0.424223,\ +184 0.304148 0.764704 0.419943,\ +185 0.311925 0.767822 0.415586,\ +186 0.319809 0.770914 0.411152,\ +187 0.327796 0.773980 0.406640,\ +188 0.335885 0.777018 0.402049,\ +189 0.344074 0.780029 0.397381,\ +190 0.352360 0.783011 0.392636,\ +191 0.360741 0.785964 0.387814,\ +192 0.369214 0.788888 0.382914,\ +193 0.377779 0.791781 0.377939,\ +194 0.386433 0.794644 0.372886,\ +195 0.395174 0.797475 0.367757,\ +196 0.404001 0.800275 0.362552,\ +197 0.412913 0.803041 0.357269,\ +198 0.421908 0.805774 0.351910,\ +199 0.430983 0.808473 0.346476,\ +200 0.440137 0.811138 0.340967,\ +201 0.449368 0.813768 0.335384,\ +202 0.458674 0.816363 0.329727,\ +203 0.468053 0.818921 0.323998,\ +204 0.477504 0.821444 0.318195,\ +205 0.487026 0.823929 0.312321,\ +206 0.496615 0.826376 0.306377,\ +207 0.506271 0.828786 0.300362,\ +208 0.515992 0.831158 0.294279,\ +209 0.525776 0.833491 0.288127,\ +210 0.535621 0.835785 0.281908,\ +211 0.545524 0.838039 0.275626,\ +212 0.555484 0.840254 0.269281,\ +213 0.565498 0.842430 0.262877,\ +214 0.575563 0.844566 0.256415,\ +215 0.585678 0.846661 0.249897,\ +216 0.595839 0.848717 0.243329,\ +217 0.606045 0.850733 0.236712,\ +218 0.616293 0.852709 0.230052,\ +219 0.626579 0.854645 0.223353,\ +220 0.636902 0.856542 0.216620,\ +221 0.647257 0.858400 0.209861,\ +222 0.657642 0.860219 0.203082,\ +223 0.668054 0.861999 0.196293,\ +224 0.678489 0.863742 0.189503,\ +225 0.688944 0.865448 0.182725,\ +226 0.699415 0.867117 0.175971,\ +227 0.709898 0.868751 0.169257,\ +228 0.720391 0.870350 0.162603,\ +229 0.730889 0.871916 0.156029,\ +230 0.741388 0.873449 0.149561,\ +231 0.751884 0.874951 0.143228,\ +232 0.762373 0.876424 0.137064,\ +233 0.772852 0.877868 0.131109,\ +234 0.783315 0.879285 0.125405,\ +235 0.793760 0.880678 0.120005,\ +236 0.804182 0.882046 0.114965,\ +237 0.814576 0.883393 0.110347,\ +238 0.824940 0.884720 0.106217,\ +239 0.835270 0.886029 0.102646,\ +240 0.845561 0.887322 0.099702,\ +241 0.855810 0.888601 0.097452,\ +242 0.866013 0.889868 0.095953,\ +243 0.876168 0.891125 0.095250,\ +244 0.886271 0.892374 0.095374,\ +245 0.896320 0.893616 0.096335,\ +246 0.906311 0.894855 0.098125,\ +247 0.916242 0.896091 0.100717,\ +248 0.926106 0.897330 0.104071,\ +249 0.935904 0.898570 0.108131,\ +250 0.945636 0.899815 0.112838,\ +251 0.955300 0.901065 0.118128,\ +252 0.964894 0.902323 0.123941,\ +253 0.974417 0.903590 0.130215,\ +254 0.983868 0.904867 0.136897,\ +255 0.993248 0.906157 0.143936) |