• API
• FAQ
• Tools
• Archive
SHARE
TWEET

# Untitled

a guest Aug 19th, 2019 65 Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. #!/usr/bin/env python2
2. # -*- coding: utf-8 -*-
3. """
4. Created on Thu Jun 27 13:05:38 2019
5.
6. @author: ltlustos
7. """
8.
9. import numpy as np
10. import pylab as pl
11. # import sympy
12. import matplotlib.pyplot as plt
13. import matplotlib.ticker
14.
15. # (x,y, x0,y0, u, X, Y, Q, Q0, s, enc, pitch, pi, sqrt2) = sympy.symbols('x y x0 y0 u X Y Q Q0 s enc pitch pi sqrt2')
16. # s_gauss_x = 1.0 / s / sqrt2 /sympy.sqrt(pi) * sympy.exp(-0.5 * (x-x0)**2 / s **2)
17. # s_gauss_y = 1.0 / s / sqrt2 / sympy.sqrt(pi) * sympy.exp(-0.5 * (y-y0)**2 / s **2)
18. # f1 = sympy.integrate(sympy.integrate(s_gauss_x*s_gauss_y , (x, -pitch/2, pitch/2)), (y, -pitch/2, pitch/2))
19.
20. datapath = '/Users/marialauraperez/Desktop/SummerCERN/Results/'
21. dataout = '/Users/marialauraperez/Desktop/SummerCERN/Results/Charge collection/'
22.
23. def pixel_charge(x0, y0, s, pitch):
24.     from scipy.special import erf
25.     from numpy import pi
26.     _sqrt2 = np.sqrt(2.0)
27.     return 0.5*pi*erf(-(1/_sqrt2)*(-pitch/2 - x0)/s)*erf(-(1/_sqrt2)*(-pitch/2 - y0)/s)/(pi*_sqrt2**2) + \
28.            0.5*pi*erf(-(1/_sqrt2)*(-pitch/2 - x0)/s)*erf((1/_sqrt2)*(pitch/2 - y0)/s)/(pi*_sqrt2**2) + \
29.            0.5*pi*erf(-(1/_sqrt2)*(-pitch/2 - y0)/s)*erf((1/_sqrt2)*(pitch/2 - x0)/s)/(pi*_sqrt2**2) + \
30.            0.5*pi*erf((1/_sqrt2)*(pitch/2 - x0)/s)*erf((1/_sqrt2)*(pitch/2 - y0)/s)/(pi*_sqrt2**2)
31.
32.
33. def normalize(arr):
34.     norm = np.linalg.norm(arr)
35.     return arr/norm
36.
37.
38. def rex_el(kv, Z, A, rhocm3):
39.     # ex_el(kv, Z, A, varargin)
40.     # um
41.     # NUCLEAR INSTRUMENTS AND METHODS IO3 (I972) 85-91;
42.     # GENERALIZED SEMIEMPIRICAL EQUATIONS
43.     # FOR THE EXTRAPOLATED RANGE OF ELECTRONS
44.     # T. TABATA, R. ITO and S. OKABE
45.     b=[.2335, 1.209, 1.78e-4, .9891, 3.01e-4, 1.468, 1.180e-2, 1.232,.109]
46.     a= [b[0]*A/Z**b[1], b[2]*Z, b[3]-b[4]*Z, b[5]-b[6]*Z, b[7]/Z**b[8]]
47.     tau= kv/511.03
48.     rex=a[0]*(np.log(1+a[1]*tau)/a[1]-a[2]*tau/(1+a[3]*tau**a[4]))/rhocm3*1e4 # um
49.     return rex
50.
51.
52. pitch_um = 55.0
53. s_um = np.array([0.71358728, 1.01030249, 1.238763, 1.43202509, 1.60287829, 1.75787494, 1.90090355, 2.03449532,
54.                  2.16040679, 2.27991584, 2.39398634, 2.50336657, 2.60865116, 2.71032194, 2.80877575, 2.90434405,
55.                  2.99730709, 3.08790435, 3.17634237, 3.26280079, 3.34743701, 3.43038988, 3.51178263, 3.59172523,
56.                  3.67031633, 3.7476448, 3.8237911, 3.89882831, 3.97282309, 4.04583643, 4.11792434, 4.18913838,
57.                  4.25952615, 4.32913174, 4.39799607, 4.46615721, 4.53365067, 4.60050961, 4.66676511, 4.73244631,
58.                  4.79758063, 4.86219387, 4.92631037, 4.98995314, 5.05314394, 5.11590341, 5.17825114, 5.24020575,
59.                  5.30178498, 5.36300571, 5.42388409, 5.48443553, 5.54467478, 5.60461599, 5.66427272, 5.72365798,
60.                  5.78278432, 5.84166378, 5.90030798, 5.95872814, 6.01693507, 6.07493925, 6.13275081, 6.19037957,
61.                  6.24783503, 6.30512645, 6.3622628, 6.41925282, 6.47610502, 6.53282769, 6.58942892, 6.64591661,
62.                  6.70229848, 6.75858209, 6.81477483, 6.87088397, 6.92691662, 6.98287977, 7.03878029, 7.09462495,
63.                  7.1504204, 7.20617322, 7.26188987, 7.31757675, 7.37324019, 7.42888644, 7.48452169, 7.54015208,
64.                  7.59578371, 7.65142261, 7.70707479, 7.76274622, 7.81844286, 7.87417062, 7.92993541, 7.98574313,
65.                  8.04159967, 8.0975109, 8.15348272])  # taken from the output of Sensor.py
66.
67. energies = np.linspace(4, 60, 113)  # spacing 0.5
68. selected_energies = np.linspace(5, 50, 10)
69. selected_weights = list()
70. all_weights = np.genfromtxt(datapath+'sigma_weights.txt')
71. for ii in range(len(energies)):
72.     if energies[ii] in selected_energies:
73.         selected_weights.append(all_weights[:, ii])
74.
75.
76. bins = np.linspace(0, 1, 100)
77. n_vals = []
78. for e in range(len(selected_energies)):
79.     n_vals.append([])
80.     for i in range(len(s_um)):
81.         X, Y = np.meshgrid(np.linspace(-pitch_um/2.0, pitch_um/2.0, 500), np.linspace(-pitch_um/2.0, pitch_um/2.0, 500))
82.         Q = pixel_charge(X, Y, s_um[i] + rex_el(e, 14, 28.08, 2.32), pitch_um)
83.         n, bins, patches = pl.hist(Q.ravel(), bins=bins, histtype='step')
84.         n_vals[e].append(n)
85. bins2 = (bins[1:] + bins[:-1]) / 2.0
86. f, ax = plt.subplots(1, 1, figsize=(8, 6))
87. for e in range(len(selected_energies)):
88.     total_spectra = np.zeros(len(n_vals[e][0]))
89.     for k in range(len(n_vals[e][0])):
90.         for j in range(len(s_um)):
91.             total_spectra[k] += n_vals[e][j][k]*selected_weights[e][j]
92.     plt.plot(bins2, normalize(total_spectra)*100, label=r'$E_\gamma$= '+str(selected_energies[e])+' keV')
93. plt.xlabel('Charge collection efficiency')
94. plt.ylabel('Percentage of occurrence')
95. plt.yscale('log')
96. plt.legend(fontsize='x-small')
97. plt.xlim(left=0.2)
98. ax.get_yaxis().set_major_formatter(matplotlib.ticker.PercentFormatter(decimals=2))
99. plt.title('Total charge collection efficiency within a pixel', fontweight='bold')
100. plt.savefig(dataout+'Charge_efficiency2.png')
RAW Paste Data
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.
Not a member of Pastebin yet?