From 6c21c214b4fdc14e45ca059d9f3c476065a2891b Mon Sep 17 00:00:00 2001 From: Gordon McCann Date: Fri, 8 Oct 2021 10:56:57 -0400 Subject: [PATCH] Added Python plotting tools using matplotlib and pickle. Updated ANASEN geometry to better reflect the actual detector setup. --- .gitignore | 1 + bin/.gitignore | 2 + bin/PyPlotViewer | 5 ++ bin/PyPlotter | 6 ++ etc/mass.txt | 3 +- include/AnasenEfficiency.h | 11 +-- input.txt | 9 +- premake5.lua | 4 +- src/Detectors/AnasenEfficiency.cpp | 4 +- src/Plotters/Python/.gitignore | 3 + src/Plotters/Python/MaskFile.py | 80 +++++++++++++++++ src/Plotters/Python/NucData.py | 70 +++++++++++++++ src/Plotters/Python/Nucleus.py | 106 ++++++++++++++++++++++ src/Plotters/Python/PyPlotter.py | 135 +++++++++++++++++++++++++++++ src/Plotters/Python/ViewPyPlots.py | 32 +++++++ 15 files changed, 455 insertions(+), 16 deletions(-) create mode 100755 bin/PyPlotViewer create mode 100755 bin/PyPlotter create mode 100644 src/Plotters/Python/.gitignore create mode 100755 src/Plotters/Python/MaskFile.py create mode 100755 src/Plotters/Python/NucData.py create mode 100755 src/Plotters/Python/Nucleus.py create mode 100755 src/Plotters/Python/PyPlotter.py create mode 100755 src/Plotters/Python/ViewPyPlots.py diff --git a/.gitignore b/.gitignore index 66caf97..cc697dc 100644 --- a/.gitignore +++ b/.gitignore @@ -12,6 +12,7 @@ .DS_Store Makefile *.make +__pycache__ ###Keep this file### diff --git a/bin/.gitignore b/bin/.gitignore index 226f3fb..888fdbf 100644 --- a/bin/.gitignore +++ b/bin/.gitignore @@ -1,3 +1,5 @@ ###keep only the directory, not the contents### * !.gitignore +!PyPlotter +!PyPlotViewer \ No newline at end of file diff --git a/bin/PyPlotViewer b/bin/PyPlotViewer new file mode 100755 index 0000000..c612b15 --- /dev/null +++ b/bin/PyPlotViewer @@ -0,0 +1,5 @@ +#!/bin/bash + +FileName=$1; + +./src/Plotters/Python/ViewPyPlots.py ${FileName} \ No newline at end of file diff --git a/bin/PyPlotter b/bin/PyPlotter new file mode 100755 index 0000000..253ddd3 --- /dev/null +++ b/bin/PyPlotter @@ -0,0 +1,6 @@ +#!/bin/bash + +InputFile=$1; +OutputFile=$2; + +./src/Plotters/Python/PyPlotter.py ${InputFile} ${OutputFile} \ No newline at end of file diff --git a/etc/mass.txt b/etc/mass.txt index da4e88f..d7fe2d5 100644 --- a/etc/mass.txt +++ b/etc/mass.txt @@ -2497,5 +2497,4 @@ 157 108 265 Hs 265 129791.799 158 108 266 Hs 266 130045.252 159 110 269 Ds 269 144751.021 - 160 110 270 Ds 270 144583.090 - + 160 110 270 Ds 270 144583.090 \ No newline at end of file diff --git a/include/AnasenEfficiency.h b/include/AnasenEfficiency.h index d59f85f..ac1da6a 100644 --- a/include/AnasenEfficiency.h +++ b/include/AnasenEfficiency.h @@ -43,14 +43,15 @@ private: const int n_qqq = 4; const double sx3_length = 0.075; const double sx3_width = 0.04; - const double barrel_gap = 0.013 + 0.049; //0.049 is base gap due to frames - const double ring1_z = sx3_length/2.0 + barrel_gap/2.0; - //const double ring2_z = -0.124 + sx3_length/2.0 + 0.0245 - barrel_gap/2.0; - const double qqq_nom_z = 0.025 + sx3_length + 0.0245 + barrel_gap/2.0; + const double barrel_gap = 0.0254; + const double sx3_frame = 0.049; //0.049 is base gap due to frames + const double ring1_z = sx3_length/2.0 + sx3_frame + barrel_gap/2.0; + const double ring2_z = (-1.0)*(barrel_gap/2.0 + sx3_length/2.0); + const double qqq_nom_z = 0.0125 + sx3_length + sx3_frame + barrel_gap/2.0; const double qqq_rinner = 0.0501; const double qqq_router = 0.0990; const double qqq_deltaphi = 1.52119; - const double qqq_z[4] = {qqq_nom_z, qqq_nom_z - 0.00828, qqq_nom_z, qqq_nom_z}; + const double qqq_z[4] = {qqq_nom_z, qqq_nom_z, qqq_nom_z, qqq_nom_z}; const double qqq_phi[4] = {5.49779, 0.785398, 2.35619, 3.92699}; const double ring_rho[12] = {0.0890601, 0.0889871, 0.0890354, 0.0890247, 0.0890354, 0.0890354, 0.0890247, 0.0890354, 0.0890354, 0.0890247, 0.0890354, 0.0890354}; const double ring_phi[12] = {0.785795, 0.262014, 6.02132, 5.49779, 4.97426, 4.45052, 3.92699, 3.40346, 2.87972, 2.35619, 1.83266, 1.30893}; diff --git a/input.txt b/input.txt index 97ca403..f3615e2 100644 --- a/input.txt +++ b/input.txt @@ -1,13 +1,14 @@ ----------Data Information---------- -OutputFile: /data1/gwm17/mask_tests/7Be12C_870keV_beam_50CD2.mask +OutputFile: /media/gordon/ANASENData/MaskData/Kinematics/7Bedp_870keV_beam_50CD2.mask SaveTree: yes SavePlots: yes ----------Reaction Information---------- -ReactionType: 1 +ReactionType: 2 Z A (order is target, projectile, ejectile, break1, break3(if pure decay is target, ejectile)) -6 12 -4 7 +1 2 4 7 +1 1 +2 4 ----------Target Information---------- Name: test_targ Layers: 1 diff --git a/premake5.lua b/premake5.lua index cd4c287..9f66766 100644 --- a/premake5.lua +++ b/premake5.lua @@ -44,8 +44,8 @@ project "RootPlot" } --User specified path to ROOT CERN libraries-- - ROOTIncludepath = "/usr/include/root" - ROOTLibpath = "/usr/lib64/root" + ROOTIncludepath = "/home/gordon/cern/root-6.22.02/root-install/include" + ROOTLibpath = "/home/gordon/cern/root-6.22.02/root-install/lib" includedirs { "include" diff --git a/src/Detectors/AnasenEfficiency.cpp b/src/Detectors/AnasenEfficiency.cpp index 560dfe7..e216578 100644 --- a/src/Detectors/AnasenEfficiency.cpp +++ b/src/Detectors/AnasenEfficiency.cpp @@ -9,7 +9,7 @@ AnasenEfficiency::AnasenEfficiency() : for(int i=0; i 0.0 else 0.0 + + px = self.vec4[0]+gfactor*bdotp*boost_vec[0]+gamma*boost_vec[0]*self.vec4[3] + py = self.vec4[1]+gfactor*bdotp*boost_vec[1]+gamma*boost_vec[1]*self.vec4[3] + pz = self.vec4[2]+gfactor*bdotp*boost_vec[2]+gamma*boost_vec[2]*self.vec4[3] + E = gamma*(self.vec4[3] + bdotp) + + self.SetVectorCart(px, py, pz, E); + + +def main(): + nuc = Nucleus(1,1) + print("First",nuc) + nuc2 = Nucleus(2,4) + print("Second", nuc2) + result = nuc + nuc2 + print("Addition", result) + result2 = nuc2 - nuc + print("Subtraction", result2) + + +if __name__ == '__main__': + main() \ No newline at end of file diff --git a/src/Plotters/Python/PyPlotter.py b/src/Plotters/Python/PyPlotter.py new file mode 100755 index 0000000..deb7042 --- /dev/null +++ b/src/Plotters/Python/PyPlotter.py @@ -0,0 +1,135 @@ +#!/usr/bin/env python3 + +import numpy as np +import matplotlib.pyplot as plt +from Nucleus import Nucleus +from MaskFile import MaskFile, MaskFileData +from NucData import Masses +import pickle +import sys + +def PlotData(inputname, outputname): + rad2deg = 180.0/np.pi + + datafile = MaskFile(inputname) + datafile.ReadHeader() + + print("MaskFile opened -- rxntype:", datafile.rxntype, "number of samples:", datafile.nsamples) + + data = MaskFileData(datafile.N_nuclei) + + ke = np.zeros((datafile.N_nuclei, datafile.nsamples)) + ke_d = np.zeros((datafile.N_nuclei, datafile.nsamples)) + theta = np.zeros((datafile.N_nuclei, datafile.nsamples)) + theta_d = np.zeros((datafile.N_nuclei, datafile.nsamples)) + phi = np.zeros((datafile.N_nuclei, datafile.nsamples)) + phi_d = np.zeros((datafile.N_nuclei, datafile.nsamples)) + detect_mask = np.ones((datafile.N_nuclei, datafile.nsamples), dtype=bool) + names = [] + for i in range(datafile.N_nuclei): + names.append(" ") + + nuc = Nucleus() + + for i in range(datafile.nsamples): + data = datafile.ReadData() + + if i == 0: + for j in range(datafile.N_nuclei): + names[j] = Masses.GetSymbol(data.Z[j], data.A[j]) + + for j in range(datafile.N_nuclei): + nuc.SetIsotope(data.Z[j], data.A[j]) + nuc.SetVectorSpher(data.theta[j], data.phi[j], data.p[j], data.E[j]) + + ke[j][i] = nuc.GetKE() + theta[j][i] = data.theta[j]*rad2deg + phi[j][i] = data.phi[j]*rad2deg + if data.dFlag[j] == True: + ke_d[j][i] = data.KE[j] + theta_d[j][i] = data.theta[j]*rad2deg + phi_d[j][i] = data.phi[j]*rad2deg + else: + detect_mask[j][i] = False + + datafile.Close() + + #Remove empty values from detection arrays + final_theta_d = theta_d[detect_mask] + final_phi_d = phi_d[detect_mask] + final_ke_d = ke_d[detect_mask] + + #figs = {} + #axes = {} + + ''' + for i in range(len(names)): + figs[i], axes[i] = plt.subplots(2,2) + ''' + fig, axes = plt.subplots(len(names)-1,4) + fig.set_size_inches(12, 12) + + for i in range(1, len(names)): + ''' + axes[i][0][0].plot(theta[i], ke[i], marker=',', linestyle='None') + axes[i][0][0].set_title(names[i]+" KE vs. Theta") + axes[i][0][0].set_xlabel(r"$\theta$ (degrees)") + axes[i][0][0].set_ylabel("KE (MeV)") + + axes[i][0][1].plot(phi[i], ke[i], marker=",", linestyle='None') + axes[i][0][1].set_title(names[i]+" KE vs. Phi") + axes[i][0][1].set_xlabel(r"$\phi$ (degrees)") + axes[i][0][1].set_ylabel("KE (MeV)") + + axes[i][1][0].plot(theta_d[i], ke_d[i], marker=',', linestyle='None') + axes[i][1][0].set_title(names[i]+" KE vs. Theta -- Detected") + axes[i][1][0].set_xlabel(r"$\theta$ (degrees)") + axes[i][1][0].set_ylabel("KE (MeV)") + + axes[i][1][1].plot(phi_d[i], ke_d[i], marker=",", linestyle='None') + axes[i][1][1].set_title(names[i]+" KE vs. Phi -- Detected") + axes[i][1][1].set_xlabel(r"$\phi$ (degrees)") + axes[i][1][1].set_ylabel("KE (MeV)") + ''' + + axes[i-1][0].plot(theta[i], ke[i], marker=',', linestyle='None') + axes[i-1][0].set_title(names[i]+" KE vs. Theta") + axes[i-1][0].set_xlabel(r"$\theta$ (degrees)") + axes[i-1][0].set_ylabel("KE (MeV)") + + axes[i-1][1].plot(phi[i], ke[i], marker=",", linestyle='None') + axes[i-1][1].set_title(names[i]+" KE vs. Phi") + axes[i-1][1].set_xlabel(r"$\phi$ (degrees)") + axes[i-1][1].set_ylabel("KE (MeV)") + + axes[i-1][2].plot(theta_d[i], ke_d[i], marker=',', linestyle='None') + axes[i-1][2].set_title(names[i]+" KE vs. Theta -- Detected") + axes[i-1][2].set_xlabel(r"$\theta$ (degrees)") + axes[i-1][2].set_ylabel("KE (MeV)") + + axes[i-1][3].plot(phi_d[i], ke_d[i], marker=",", linestyle='None') + axes[i-1][3].set_title(names[i]+" KE vs. Phi -- Detected") + axes[i-1][3].set_xlabel(r"$\phi$ (degrees)") + axes[i-1][3].set_ylabel("KE (MeV)") + + plt.tight_layout() + plt.show(block=True) + + print("Writing figure to file:", outputname) + with open(outputname, "wb") as outfile: + pickle.dump(fig, outfile) + outfile.close() + + print("Finished.") + +def main(): + if len(sys.argv) == 3: + PlotData(sys.argv[1], sys.argv[2]) + else: + print("Unable to run PyPlotter, incorrect number of arguments! Need an input datafile name, and an output plot file name") + +if __name__ == '__main__': + main() + + + diff --git a/src/Plotters/Python/ViewPyPlots.py b/src/Plotters/Python/ViewPyPlots.py new file mode 100755 index 0000000..3ed324d --- /dev/null +++ b/src/Plotters/Python/ViewPyPlots.py @@ -0,0 +1,32 @@ +#!/usr/bin/env python3 + +import matplotlib.pyplot as plt +import numpy as np +import pickle +import sys + +def SetManager(figure): + dummy = plt.figure() + manager = dummy.canvas.manager + manager.canvas.figure = figure + figure.set_canvas(manager.canvas) + +def ViewPyPlots(filename): + + figure = pickle.load(open(filename, "rb")) + SetManager(figure) + + figure.set_size_inches(12, 12) + plt.tight_layout() + + plt.show(block=True) + + +def main(): + if len(sys.argv) == 2: + ViewPyPlots(sys.argv[1]) + else: + print("Unable to run ViewPyPlots, incorrect number of commandline arguments -- requires an input pickle file") + +if __name__ == '__main__': + main() \ No newline at end of file