From 0397ef71f88a424fb23076f936e20142c9818c6a Mon Sep 17 00:00:00 2001 From: liuyihui Date: Thu, 1 Jun 2023 15:15:10 +0800 Subject: [PATCH] feat: data process script --- .gitignore | 1 + init_vis.mac | 8 +-- integral.cpp | 69 ++++++++++++++++++------- run.mac | 2 +- script/draw | Bin 0 -> 17864 bytes script/draw.cpp | 76 +++++++++++++++++++++++++++ script/draw.py | 127 ++++++++++++++++++++++++++++++++++++++++++++++ script/eff.py | 32 ++++++++++++ script/eff.txt | 80 +++++++++++++++++++++++++++++ script/eff2.py | 24 +++++++++ script/eff2.txt | 15 ++++++ script/eff3.py | 24 +++++++++ script/ej-200.txt | 21 ++++++++ script/inner.txt | 15 ++++++ script/main.py | 16 ++++++ script/outer.txt | 15 ++++++ script/roi.py | 27 ++++++++++ script/roi.txt | 80 +++++++++++++++++++++++++++++ script/roi2.py | 28 ++++++++++ script/roi2.txt | 15 ++++++ 20 files changed, 653 insertions(+), 22 deletions(-) create mode 100755 script/draw create mode 100644 script/draw.cpp create mode 100644 script/draw.py create mode 100644 script/eff.py create mode 100644 script/eff.txt create mode 100644 script/eff2.py create mode 100644 script/eff2.txt create mode 100644 script/eff3.py create mode 100644 script/ej-200.txt create mode 100644 script/inner.txt create mode 100644 script/main.py create mode 100644 script/outer.txt create mode 100644 script/roi.py create mode 100644 script/roi.txt create mode 100644 script/roi2.py create mode 100644 script/roi2.txt diff --git a/.gitignore b/.gitignore index d58425e..be48023 100644 --- a/.gitignore +++ b/.gitignore @@ -15,3 +15,4 @@ test/ # Data *.root *.data +data/ diff --git a/init_vis.mac b/init_vis.mac index 2ed5342..b7f97a6 100644 --- a/init_vis.mac +++ b/init_vis.mac @@ -1,7 +1,9 @@ /control/verbose 1 -/run/verbose 2 -/event/verbose 0 -/tracking/verbose 0 +/run/verbose 1 +/event/verbose 1 +/tracking/verbose 1 + +/Gene/PrimaryGA/SetBeamEnergy 1 MeV # Initialize kernel /run/initialize diff --git a/integral.cpp b/integral.cpp index 15380fc..5008b91 100644 --- a/integral.cpp +++ b/integral.cpp @@ -1,39 +1,72 @@ #include +#include #include #include -void integral(TString file) { +int integral(TString file) { double e[24]; int id[24], N[1]; - int sum1, sum2; + int sum1, sum2, sum3; TFile f(file); - if (f.IsZombie()) { - printf("cannot open root file\n"); - return; - } TTree* t = (TTree*)f.Get("SimData1"); - if (!t) { - printf("cannot find tree!\n"); - return; - } t->SetBranchAddress("detE", e); t->SetBranchAddress("detId", id); t->SetBranchAddress("Ndets", N); int n = t->GetEntries(); // printf("n=%d\n", n); int i, j; - sum1 = sum2 = 0; + sum1 = sum2 = sum3 = 0; for (i = 0; i < n; i++) { t->GetEntry(i); for (j = 0; j < N[0]; j++) { - if (e[j] > 0.18 && id[j] < 12) { - sum1++; - } - if (e[j] > 0.18 && id[j] > 11) { - sum2++; - } + if (e[j] > 0.18 && id[j] < 12) sum1++; + if (e[j] > 0.18 && id[j] > 11) sum2++; + if (e[j] >= 0.64 && e[j] <= 0.8) sum3++; } } - printf("%d,%d,%d,%f\n", sum1 + sum2, sum1, sum2, 1.0 * sum1 / sum2); + // printf("%d, %d, %d, %d, %f\n", n, sum1 + sum2, sum1, sum2, 1.0 * sum1 / sum2); + // return sum1 + sum2; + return sum3; +} + +void process() { + double or1[5] = {10, 11, 12, 13, 14}; + double or2[5] = {12.14, 13.35, 14.57, 15.79, 17}; + double or3[5] = {14.28, 15.71, 17.14, 18.57, 20}; + // int energy[8] = {600, 650, 700, 750, 800, 850, 900, 950}; + + for (int i = 0; i < 5; i++) { + TString file = + TString::Format("data/7/%d/Simulation_RunNo_650kev_70mm_%dmm.root", int(or1[i]), int(or1[i] * 10)); + printf("%d\n", integral(file)); + } + + for (int i = 0; i < 5; i++) { + TString file = + TString::Format("data/8.5/%.2f/Simulation_RunNo_650kev_85mm_%dmm.root", or2[i], int(or2[i] * 10)); + printf("%d\n", integral(file)); + } + + for (int i = 0; i < 5; i++) { + TString file = + TString::Format("data/10/%.2f/Simulation_RunNo_650kev_100mm_%dmm.root", or3[i], int(or3[i] * 10)); + printf("%d\n", integral(file)); + } + + // for (int i = 0; i < 5; i++) { + // for (int j = 0; j < 8; j++) { + // TString file = TString::Format("data/7/%d/Simulation_RunNo_%dkev_70mm_%dmm.root", int(or1[i]), energy[j], + // int(or1[i] * 10)); + // printf("%d\n", integral(file)); + // } + // } + + // for (int i = 0; i < 5; i++) { + // for (int j = 0; j < 8; j++) { + // TString file = TString::Format("data/8.5/%.2f/Simulation_RunNo_%dkev_85mm_%dmm.root", or2[i], energy[j], + // int(or2[i] * 10)); + // printf("%d\n", integral(file)); + // } + // } } diff --git a/run.mac b/run.mac index f238e54..dc6d816 100644 --- a/run.mac +++ b/run.mac @@ -6,4 +6,4 @@ /run/initialize /gun/particle neutron -/control/loop loop.mac energy 0.1 5 0.1 +/control/loop loop.mac energy 0.5 1 0.05 diff --git a/script/draw b/script/draw new file mode 100755 index 0000000000000000000000000000000000000000..7b999d4d1a21955151e4ebfa8912ae2f68c0e24a GIT binary patch literal 17864 zcmeHP4RBP~b-vQFAbzx35o9o~%{n0sm}rqeVBrLNS7>47NLUe&NDlB=?Jj6N+Fh~x z77MpH^2qj}I%;K4JtiJ|DtDTO$UC@S&2Lfj(ef}bNXB@ZY9sVe8n3u%SMU7%#Qj4DOw?^ikw zr?zA2)pkrdpX7!Vu4mPksyGwCmJ%PTH*XnNk~6cD!lEuBW|rZ8zYPZiJJ@ z6fBr>d!L3K<>hxRJj@-sz1hZLT<4Q%nQmZhZ=!o+!`j|hZEqr#>#yx^+*rG@!Jke0 z*GUbEi`t+$wQc(@0XO4>=`_l*SfY87KRbBiSN47R?9ZRPaL3c(M@Jv}i^@i&iy+R|4Fju*rF{`xyaWB7%nb?8^DU;AD>14W3%8}X2DO* z;{SzN@XNE{D}ej(ahCd7+BIguDPKOP6Y!MvVtM!ka34O-@-%>A`t4}@;@P8KE5;(KMWxjfr&DfJVHbC6zFn>qCbmPR`r6cI-&+Fha78-fYY7 z=0Ca<-ZvQSpNbZ1V2f%1F-OxmQy@18UsJ`TU6F>4ZEz zb=KnPLSlOAq{Y*Td3x#_7EdSS>8a-|o=(KmQ;%Cboq(sO9{)mImIryZ5kK6hBw_1kI?hOz9OL*|#F0}7# z9d`dddK4b6d>$PQ4{x4A|MN8+F!}We*u$uVhus+>n=Y9(7#cs-F*Nh@6ICk)(0g%8 zQ|Pxy{S8e-R(}v4x)473_FdtVmt5hpSHfpMFqgoATRW)8Pwba@VtX3Dfz2TduE?#~ z6&~FD0zG&L54~wF3O}~_-@#0N@q>JRGKOTo;=UPrWqV+2wNLg__d`gPU9by3bbbUL zrmLzS=?6Sr<#mtaf(Lx)D3~_-a=$2%u9Jk0Sm<#ZddxzfwxOeB{J3^^91K}kyPrbu zh=EbH`w7iJx7tnPN)qT+yN@cy-J|HRimy{fF$OIe@cyb*r>CpB@Ui$a(^Z}Lc-?EI zVb>49OQRFv$3iFehpSTKlYe?OpGWP%@$gV+oOBneM%v(5DU(MM?G2$;?RF`7<>gu# zq5GnwvR(JLa0NYGwfp4co$v~?%D&G$p=da}R;h`U33jO&S znisExhyG#m7XYUgkLS%Y+edi!;4_$8a+7afKPhv6lp^MXp} zi9ayNG>VKGd4*lB(I*b4urw;`zX1#Nn${mqFJk2lDQw7nL}e*?1dOTnm`WU zfQR^Tvbl`(Q-Uu_D!X3uEsZ^|FxlnF1=s=?&4fLS1jZEBQ8YIg} z@-d;rGHHa+G1GZ+a^g2LU7=Zm(VjW@Ywp#}btrEtxnG9~#j;s~>gb_)qGkgQn)9o0 z(7XpNbVxJ5fFlLd!;7Qp$2gbGB=(1)D_A?l;E^lQ&##)Ux~=UB_2)En^hj2| zEbDtuPCkhvd@yvS?3b^fAd^DJLr1O*r>-@m>v@e_P&6$aG>Q<7!_FU-2i#H7i-aC*?ex=+S zm~X7R@3pHN0*5~H>faex8UxFYxD$`=9u6p9{M27ae&At$Q2G7MUF=`!bp^_9a(yBO zW+6VKzsl#2fd($*^W&iHpyxm*KnFlyeJ7v40Q&l`^LhFm?VbP3=j%Z4#$^Wm{`M4T zALtn92&fO&GUq_4&~LF+${y$tW&Pf=n-REPo)x!HL2E<*<@3?LK$`2Dv_Pg=vho1d+6~aw?`he5@Fc)g;3;3Kxd-LhJ zNXX-Tw0ukDT$gYmBpJ)s0DtFVKEHy7#**Jq&hj4v{`R~1d>5!CFK;KYv_B7c8uGRR zd9x*7gnq9;?ryZpmAPmJM!x|sgFH_1tsd{6l(%@6{jn?L@eR)jc~*VF-R!A*bZ*$w zI5=;cCy??q20e8_&#ESmugSBl$>VMERLI|$O!x`F4?W%~l~M*u87O6-lz~zPN*O3+ zpp=192L3-}fZx;N_pUsw2@r_)-#K&go#9UtAFQPJUu8pT@Q z+D!E5_ni2ylJ4KA(ESh9?al&ot zEpCkUTs2{gAp9}Tp3ZBP>_&mY$y%4VuJP!)-`GV(9A-+KH zJT1f*TJyOO$I)G6oLtyN;VX)J#CTD>N}MT*d#(K%h3!@gUOx(Pyj<$DBdYEaAF|%Z zw@Y0(dMAtO*I4UOq5cil`>BQaQc>(UaA8&9b)!%pt6sp4sJhGcKDAxy61*zTptO9< z5%1tvtinR;Z8-*I#Qu^gFptK^9r&#p=k*<|EVtvk4DqxF6;JI^*@{n@xHb<-;M6Xk zH%OPtVSLK4s`C0SL3|yYN&Us5!iNe9xBJ(!zFdqs;(1Nti-p&R3JR~AucJ=! z{75gz@|h7~a7Z|x*ZBNH?G%gW1_;V9PkQB3F(}zeG5wX~sK|VyKKlMN%v#fMayUK|PLL%1pc|an0CB zZ8_pEzfthI&!UGB=i>5f*Nr=%Uo1bZPxR=^R{QLBAw^j%A;4J#5 zX2IV8PW^I@*IR@m&SQ=^FKRpxPz`#KFDeTV4~@GBNY?slh1?52C$t|HES37puW9!q zvR|wV*}QAfLfA@%K^RsEZ@;C~O?cdfek*ID%6pGAKO4uoR!O++(h*397z zfvA|-A7hwFBT74Dve@|%OB)Az)7_C?BW9*EStF9`$L^72UvJ!u$NY_T8|#ZJ(e9Ci z5y@mChmCm3%p4Z`Gm&K6h~<*W!!WUN1~knA)mSESzdzd7CwIlP)^F-d#G~>16WO>C zO{cKOCKtsX6>Q3BZ8VHTn#{qruIo#HG}IfkQ^mw~912!%BmZtf69bJzA@Tgw*63PQ(D zAv++qA6YHjvvuo^5He-Ame;-sgk4kWCvs$Rf5)(JrFG zw)z?tr_my3cQRsVrBSgK^>QQ^l1H?T72{!CgPbN+abUXFTV z@W!0-q1Gj3!F4nhNhaV-nPajDKi0$~*0~~!^>sjWnQvZAJZoyBmA!=CRqnfh2!Q>XoBf#avn zG{?`+S(vWUcKn>h8Gi(Fdd9-`eE-UHjFO8&eWWVe^Yf4dFzPbf^L;K;evShbS+qN5 z#X;2Q8426-{V-F`AIER$X%7#A(lZyv`M#OyKCMUbQ=A;XhQ%-{lvB3n`)j5?9Wb}g z_T2xcw0)Cy$j_CSp3^wlla15(6FC=CJ4I1WZTRP*e8D^xF>mefm0_X@~~SR@jcA@1tVRA3t|s%Fi=d->)YA>&mghe^Yad-M@eWa&i3DgjGi6n zWyiXY=l?t4%op1;eFwHxHhX@K+g=M0h4mCS!!v#r0-C$rK0jyS|0BWwFM@}W?KrNl z!j`(o_WV4iN8gAYgUnvop6Of|+wDa_Yu72t?XVBFW1ikgA&lD(1eC<9?fc2Otpsd! zs?z-wsc3ae%Jyq4L{-?Go$}nVs|K{gJ{|DQ5U$6EzjxJ`};$L+E% S9d5= 1.0e6) + thetalabr = pi / 2.0; + else + thetalabr = std::atan(tg_thetalabr); + if (thetalabr < 0.0) thetalabr = pi + thetalabr; + + double p4_cmx = p4_cm * std::sin(thetacmr); + double p4_cmz = p4_cm * std::cos(thetacmr); + double p4_labx = p4_cmx; + double p4_labz = vgamma * (p4_cmz + vbeta * (t4_cm + mass[3])); + double p4_lab = std::sqrt(p4_labx * p4_labx + p4_labz * p4_labz); + double e4_lab = sqrt(p4_lab * p4_lab + mass[3] * mass[3]) - mass[3]; + + // std::cout << thetalabr << ", " << e4_lab << std::endl; + std::cout << e4_lab << std::endl; +} + +int main() { + SetBeamEnergy(1); + for (int i = 0; i < 100; i++) GeneratePrimaries(i * pi / 100.0); + return 0; +} diff --git a/script/draw.py b/script/draw.py new file mode 100644 index 0000000..7526c8e --- /dev/null +++ b/script/draw.py @@ -0,0 +1,127 @@ +import numpy as np +from matplotlib import pyplot as plt + +mu = 931.4940954 # MeV/c^2 + + +class Particle(object): + """Define a Particle + + Parameters + ---------- + A : int + mass number + m : float + mass + """ + + def __init__(self, A, m) -> None: + self.A = A + self.m = m + + +class Reaction(object): + """A Nuclear Reaction A(a,b)B + target nuclei (A) is assumed to be quiescent + + Parameters + ---------- + a/A/b/B : Particle + Particle a, A, b, B + """ + + def __init__(self, a, A, b, B) -> None: + self.a: Particle = a + self.A: Particle = A + self.b: Particle = b + self.B: Particle = B + + @property + def Q(self): + """Returns Q value of the reaction""" + return (self.A.m + self.a.m - self.B.m - self.b.m) * mu + + @property + def threshold(self): + """Returns threshold energy of the reaction""" + if self.Q < 0: + return -(self.A.m + self.a.m) / self.A.m * self.Q + return 0 + + def eb(self, theta, ke): + """Returns energy of b in the lab frame + + Parameters + ---------- + theta : float + angle between a and b in the lab frame in rad + ke : float + kinetic energy of a in MeV + """ + i1 = np.sqrt(self.a.A * self.b.A * ke) / (self.B.A + self.b.A) * np.cos(theta) + i2 = (self.B.A - self.a.A) / (self.B.A + self.b.A) + (self.a.A * self.b.A) / ((self.B.A + self.b.A) ** 2) * (np.cos(theta) ** 2) + i3 = self.B.A * self.Q / (self.B.A + self.b.A) + return (i1 + np.sqrt(i2 * ke + i3)) ** 2 + + def eB(self, theta, ke): + """Returns energy of B in the lab frame + + Parameters + ---------- + theta : float + angle between a and b in the lab frame in rad + ke : float + kinetic energy of a in MeV + """ + eb = self.eb(theta, ke) + i1 = self.a.A / self.B.A * ke + self.b.A / self.B.A * eb + i2 = 2 * np.sqrt(self.a.A * self.b.A * ke * eb) / self.B.A + return i1 - i2 * np.cos(theta) + + +neutron = Particle(1, 1.008665) +proton = Particle(1, 1.007825) +tritium = Particle(3, 3.016049) +he_3 = Particle(3, 3.016029) +he_4 = Particle(4, 4.002603) +ne_22 = Particle(22, 21.991385) +mg_25 = Particle(25, 24.985837) + +theta = np.linspace(0, np.pi, 100) + +r1 = Reaction(he_4, ne_22, neutron, mg_25) +res = r1.eb(theta, 1) + +plt.plot(theta, res) +plt.title(r"$^{22}$Ne($\alpha$,n)$^{25}$Mg, $T_\alpha$=0.8 MeV") +plt.ylabel(r"$T_n$ (MeV)") +plt.xlabel(r"$\theta$ (rad)") +plt.xticks([0, np.pi / 4, np.pi / 2, np.pi * 3 / 4, np.pi], [r"$0$", r"$\pi/4$", r"$\pi/2$", r"$3\pi/4$", r"$\pi$"]) +plt.xlim(0, np.pi) +plt.show() + +exit(0) + +energy = np.linspace(np.min(res), np.max(res), 100) + +r2 = Reaction(he_3, neutron, proton, tritium) + +# res = r2.eb(theta, 0.0114) +# res2 = r2.eB(theta, 0.0114) + +res = r2.eb(np.pi / 4, energy) +res2 = r2.eB(np.pi / 4, energy) + +fig, ax = plt.subplots() + +ax.plot(energy, res, color="C0") +ax.set_title(r"$^3$He(n,p)$^3$H, $\theta$=$\pi/4$") +ax.set_xlabel(r"$T_n$ (MeV)") +ax.set_ylabel(r"$T_p$ (MeV)") + +ax1 = ax.twinx() +ax1.plot(energy, res2, color="C1") +ax1.set_ylabel(r"$T_t$ (MeV)") + +plt.tight_layout() +plt.show() diff --git a/script/eff.py b/script/eff.py new file mode 100644 index 0000000..edff3da --- /dev/null +++ b/script/eff.py @@ -0,0 +1,32 @@ +import numpy as np +from matplotlib import pyplot as plt + +or1 = [10, 11, 12, 13, 14] +or2 = [12.14, 13.35, 14.57, 15.79, 17] +energy = [600, 650, 700, 750, 800, 850, 900, 950] +markers = ["^", "+", "s", "D", "."] + +data = np.loadtxt("test/eff.txt") +data = data.reshape(10, 8) +data = data / 10000 + +plt.figure(dpi=400) + +ax = plt.subplot(1, 1, 1) +for k in range(5): + ax.scatter(energy, data[k], label="{:d}cm".format(or1[k]), marker=markers[k]) +ax.legend(loc="lower left") +ax.set_xlabel("Energy [keV]", fontsize=15) +ax.set_ylabel("Efficiency", fontsize=15) + +# ax = plt.subplot(2, 1, 2) +# for k in range(5): +# ax.scatter(energy, data[k + 5], label="{:.2f}cm".format(or2[k]), marker=markers[k]) +# ax.legend(loc="lower left") +# ax.set_xlabel("Energy [keV]", fontsize=15) +# ax.set_ylabel("Efficiency", fontsize=15) + +plt.tight_layout() +plt.subplots_adjust() +# plt.show() +plt.savefig("/home/fox/Documents/Thesis/figures/chap3/Fig3.png") diff --git a/script/eff.txt b/script/eff.txt new file mode 100644 index 0000000..9fa32b8 --- /dev/null +++ b/script/eff.txt @@ -0,0 +1,80 @@ +5389 +5171 +4988 +4774 +4722 +4599 +4533 +4401 +5726 +5605 +5383 +5255 +5089 +5022 +4959 +4824 +5880 +5750 +5613 +5496 +5341 +5233 +5190 +5094 +5885 +5789 +5669 +5515 +5480 +5345 +5317 +5121 +5751 +5649 +5455 +5535 +5436 +5266 +5141 +5253 +5766 +5669 +5512 +5378 +5405 +5283 +5127 +5085 +5819 +5720 +5632 +5544 +5472 +5485 +5456 +5162 +5650 +5642 +5507 +5456 +5437 +5366 +5368 +5361 +5567 +5492 +5435 +5347 +5350 +5248 +5275 +5182 +5412 +5353 +5310 +5181 +5219 +5121 +5061 +5055 \ No newline at end of file diff --git a/script/eff2.py b/script/eff2.py new file mode 100644 index 0000000..498a1a7 --- /dev/null +++ b/script/eff2.py @@ -0,0 +1,24 @@ +import numpy as np +from matplotlib import pyplot as plt + +or1 = [10, 11, 12, 13, 14] +markers = ["^", "+", "s", "D", "."] + +data = np.loadtxt("test/eff2.txt") +data = data.reshape(3, 5) +data = data / 10000 + +fig = plt.figure(dpi=400) + +ax = plt.subplot(1, 1, 1) +for k in range(5): + ax.scatter([7, 8.5, 10], data[:, k], label=r"$R_o/R_i$ = {:.2f}".format(or1[k] / 7), marker=markers[k]) +ax.legend(loc="lower left") +ax.set_xlabel(r"$R_i$ [cm]", fontsize=15) +ax.set_ylabel("Efficiency", fontsize=15) +ax.set_xticks([7, 8.5, 10]) + +plt.tight_layout() +plt.subplots_adjust() +# plt.show() +plt.savefig("/home/fox/Documents/Thesis/figures/chap3/Fig5.png") diff --git a/script/eff2.txt b/script/eff2.txt new file mode 100644 index 0000000..7909d6d --- /dev/null +++ b/script/eff2.txt @@ -0,0 +1,15 @@ +5171 +5605 +5750 +5789 +5649 +5669 +5720 +5642 +5492 +5353 +5165 +5030 +4957 +4831 +4642 \ No newline at end of file diff --git a/script/eff3.py b/script/eff3.py new file mode 100644 index 0000000..2d4dfdc --- /dev/null +++ b/script/eff3.py @@ -0,0 +1,24 @@ +import numpy as np +from matplotlib import pyplot as plt + +or1 = [10, 11, 12, 13, 14] +markers = ["^", "+", "s", "D", "."] + +data = np.loadtxt("test/inner.txt") +data = data.reshape(3, 5) +data = data / 10000 + +fig = plt.figure(dpi=400) + +ax = plt.subplot(1, 1, 1) +for k in range(5): + ax.scatter([7, 8.5, 10], data[:, k], label=r"$R_o/R_i$ = {:.2f}".format(or1[k] / 7), marker=markers[k]) +ax.legend(loc="lower right") +ax.set_xlabel(r"$R_i$ [cm]", fontsize=15) +ax.set_ylabel("Efficiency", fontsize=15) +ax.set_xticks([7, 8.5, 10]) + +plt.tight_layout() +plt.subplots_adjust() +# plt.show() +plt.savefig("/home/fox/Documents/Thesis/figures/chap3/Fig6(a).png") diff --git a/script/ej-200.txt b/script/ej-200.txt new file mode 100644 index 0000000..c5188ae --- /dev/null +++ b/script/ej-200.txt @@ -0,0 +1,21 @@ +395.00000, 0.00000 +400.00000, 0.00926 +403.19005, 0.03801 +408.51147, 0.19876 +410.83959, 0.30444 +412.75197, 0.40119 +418.82171, 0.79712 +421.48241, 0.93703 +422.81277, 0.97870 +425.00000, 1.00000 +427.63530, 0.98168 +433.37245, 0.89833 +440.10736, 0.74651 +446.75913, 0.60064 +460.06266, 0.43394 +465.71666, 0.37440 +472.45158, 0.26426 +480.10111, 0.17197 +490.16190, 0.10053 +499.97326, 0.05885 +515.00000, 0.00000 \ No newline at end of file diff --git a/script/inner.txt b/script/inner.txt new file mode 100644 index 0000000..df7a5aa --- /dev/null +++ b/script/inner.txt @@ -0,0 +1,15 @@ +2295 +2811 +3105 +3346 +3515 +3543 +3887 +4070 +4253 +4331 +3835 +3960 +4153 +4236 +4213 \ No newline at end of file diff --git a/script/main.py b/script/main.py new file mode 100644 index 0000000..e3eb57a --- /dev/null +++ b/script/main.py @@ -0,0 +1,16 @@ +import numpy as np + +data = np.loadtxt("ej-200.txt", delimiter=",") + +c = 299792458 +h = 6.62606957 * 1e-34 +hc = h * c / (1.6 * 1e-19) + +L = len(data) +for k in range(L): + print(str(round(hc * 1e9 / data[L - k - 1][0], 3)) + " * eV", end=", ") +print() + +for k in range(L): + print("{:.5f}".format(round(data[L - k - 1][1], 5)), end=", ") +print() diff --git a/script/outer.txt b/script/outer.txt new file mode 100644 index 0000000..e05c2ca --- /dev/null +++ b/script/outer.txt @@ -0,0 +1,15 @@ +2876 +2794 +2645 +2443 +2134 +2126 +1833 +1572 +1239 +1022 +1330 +1070 +804 +595 +429 \ No newline at end of file diff --git a/script/roi.py b/script/roi.py new file mode 100644 index 0000000..8e7619c --- /dev/null +++ b/script/roi.py @@ -0,0 +1,27 @@ +import numpy as np +from matplotlib import pyplot as plt + +or1 = [10, 11, 12, 13, 14] +or2 = [12.14, 13.35, 14.57, 15.79, 17] +energy = [600, 650, 700, 750, 800, 850, 900, 950] +markers = ["^", "+", "s", "D", "."] + +data = np.loadtxt("test/eff.txt") +data = data.reshape(10, 8) + +roi = np.loadtxt("test/roi.txt") +roi = roi.reshape(10, 8) +roi = 1 - roi / data + +plt.figure(dpi=400) + +for k in range(5): + plt.scatter(energy, roi[k], label="{:d}cm".format(or1[k]), marker=markers[k]) +plt.xlabel("Energy [keV]", fontsize=15) +plt.ylabel(r"Energy Leakage Rate $\epsilon$", fontsize=15) +plt.legend(loc="lower left") + +plt.tight_layout() +plt.subplots_adjust() +# plt.show() +plt.savefig("/home/fox/Documents/Thesis/figures/chap3/Fig7(a).png") diff --git a/script/roi.txt b/script/roi.txt new file mode 100644 index 0000000..0812ed9 --- /dev/null +++ b/script/roi.txt @@ -0,0 +1,80 @@ +3973 +3798 +3656 +3505 +3442 +3391 +3363 +3277 +4219 +4083 +3952 +3875 +3694 +3634 +3680 +3532 +4233 +4283 +4186 +4057 +3921 +3818 +3758 +3745 +4342 +4214 +4178 +4024 +3992 +3954 +3888 +3721 +4202 +4161 +4027 +4012 +3971 +3855 +3760 +3819 +4252 +4194 +4073 +3957 +3945 +3873 +3761 +3638 +4296 +4125 +4112 +4117 +4037 +3986 +3951 +3748 +4144 +4138 +4041 +4001 +3965 +3924 +3951 +3880 +4067 +3990 +3938 +3902 +3933 +3812 +3806 +3772 +3922 +3903 +3913 +3755 +3813 +3772 +3612 +3724 \ No newline at end of file diff --git a/script/roi2.py b/script/roi2.py new file mode 100644 index 0000000..f86db10 --- /dev/null +++ b/script/roi2.py @@ -0,0 +1,28 @@ +import numpy as np +from matplotlib import pyplot as plt + +or1 = [10, 11, 12, 13, 14] +or2 = [12.14, 13.35, 14.57, 15.79, 17] +energy = [600, 650, 700, 750, 800, 850, 900, 950] +markers = ["^", "+", "s", "D", "."] + +data = np.loadtxt("test/eff2.txt") +data = data.reshape(3, 5) + +roi = np.loadtxt("test/roi2.txt") +roi = roi.reshape(3, 5) +roi = 1 - roi / data + +plt.figure(dpi=400) + +for k in range(5): + plt.scatter([7, 8.5, 10], roi[:, k], label=r"$R_o/R_i$ = {:.2f}".format(or1[k] / 7), marker=markers[k]) +plt.xlabel(r"$R_i$ [cm]", fontsize=15) +plt.xticks([7, 8.5, 10]) +plt.ylabel(r"Energy Leakage Rate $\epsilon$", fontsize=15) +plt.legend(loc="lower right") + +plt.tight_layout() +plt.subplots_adjust() +# plt.show() +plt.savefig("/home/fox/Documents/Thesis/figures/chap3/Fig7(b).png") diff --git a/script/roi2.txt b/script/roi2.txt new file mode 100644 index 0000000..d466fe3 --- /dev/null +++ b/script/roi2.txt @@ -0,0 +1,15 @@ +3798 +4083 +4283 +4214 +4161 +4194 +4125 +4138 +3990 +3903 +3745 +3638 +3624 +3532 +3410 \ No newline at end of file