Skip to content
Snippets Groups Projects
Commit a93bf3b1 authored by Tooba Ali's avatar Tooba Ali
Browse files

include Sigma and eSigma and Kinematic Coverage script in... 

include Sigma and eSigma and Kinematic Coverage script in benchmarks/dis/analysis/kinematics_correlations.py
parent bd1a4390
No related branches found
No related tags found
1 merge request!202construct kinematic coverage plot and include Sigma and eSigma correlation in kinematics_correlation.py
This commit is part of merge request !202. Comments created here will be created in the context of that merge request.
Hide whitespace changes
Inline Side-by-side
benchmarks/dis/analysis/kinematics_correlations.py
+ 69
− 17
View file @ a93bf3b1
...@@ -14,18 +14,20 @@ import argparse ...@@ -14,18 +14,20 @@ import argparse
parser = argparse.ArgumentParser() parser = argparse.ArgumentParser()
parser.add_argument('--rec_file', type=str, help='Reconstructed track file.') parser.add_argument('--rec_file', type=str, help='Reconstructed track file.')
parser.add_argument('--ebeam', type=float, help='Electron beam energy.') parser.add_argument('--config', type=str, help='Momentum configuration.')
parser.add_argument('--pbeam', type=float, help='Proton (or ion) beam energy.') parser.add_argument('--nevents', type=float, help='Number of events to process.')
parser.add_argument('--minq2', type=float, help='Minimum four-momentum transfer squared Q2.') parser.add_argument('--results_path', type=str, help='Output directory.')
parser.add_argument('-o', dest='outdir', default='results/dis/', help='Output directory.')
args = parser.parse_args() args = parser.parse_args()
kwargs = vars(args) kwargs = vars(args)
rec_file = args.rec_file rec_file = args.rec_file
minq2 = int(args.minq2) Nevents = int(args.nevents)
k = int(args.ebeam) r_path = args.results_path #Path for output figures and file.
p = int(args.pbeam) Dconfig = 'epic' + args.config.split('_epic')[1].strip() #Detector config
config = args.config.split('_epic')[0].strip()
minq2 = int(config.split('=')[1].strip()[0])
k = int(config.split('x')[0].strip()) #ebeam Electron beam energy
p = int(config.split('x')[1].split('_')[0].strip()) #pbeam Proton (or ion) beam energy
#logarithmically spaced bins in Q2 and 𝑥 #logarithmically spaced bins in Q2 and 𝑥
...@@ -33,8 +35,10 @@ Q2_bins = [0.40938507,0.64786184,1.025353587,1.626191868,2.581181894,4.091148983 ...@@ -33,8 +35,10 @@ Q2_bins = [0.40938507,0.64786184,1.025353587,1.626191868,2.581181894,4.091148983
x_bins = [4.09385E-05,6.47862E-05,0.000102535,0.000162619,0.000258118,0.000409115,0.000647862,0.001025354,0.001626192,0.002581182,0.004091149,0.006478619,0.010253536,0.016261919,0.025811819,0.04091149,0.064786187,0.10253536,0.162619187,0.25811819,0.409114898,0.647861868,1.025257131] x_bins = [4.09385E-05,6.47862E-05,0.000102535,0.000162619,0.000258118,0.000409115,0.000647862,0.001025354,0.001626192,0.002581182,0.004091149,0.006478619,0.010253536,0.016261919,0.025811819,0.04091149,0.064786187,0.10253536,0.162619187,0.25811819,0.409114898,0.647861868,1.025257131]
##########################################################################################
#function to construct Q2 correlation plots #function to construct Q2 correlation plots
##########################################################################################
def Q2correlation(minq2,method): #minq2 can be 1,10,100, or 1000; method can be 'e','DA', or 'JB' def Q2correlation(minq2,method): #minq2 can be 1,10,100, or 1000; method can be 'e','DA', or 'JB'
Q2values_Y = method_Q2values_dict['{}'.format(method)] #Q2 values of the given method, that are mapped onto the y axis Q2values_Y = method_Q2values_dict['{}'.format(method)] #Q2 values of the given method, that are mapped onto the y axis
...@@ -85,12 +89,13 @@ def Q2correlation(minq2,method): #minq2 can be 1,10,100, or 1000; method can be ...@@ -85,12 +89,13 @@ def Q2correlation(minq2,method): #minq2 can be 1,10,100, or 1000; method can be
plt.xlabel('$Q^2$ [$GeV^2$] Truth') plt.xlabel('$Q^2$ [$GeV^2$] Truth')
plt.ylabel('$Q^2$ [$GeV^2$] {}'.format(method_dict['{}'.format(method)])) plt.ylabel('$Q^2$ [$GeV^2$] {}'.format(method_dict['{}'.format(method)]))
plt.title('{} $Q^2$ correlation {}x{} $minQ^2=${}$GeV^2$'.format(method_dict['{}'.format(method)],k,p,minq2)) plt.title('{} $Q^2$ correlation {}x{} $minQ^2=${}$GeV^2$'.format(method_dict['{}'.format(method)],k,p,minq2))
plt.show() plt.savefig(os.path.join(r_path, 'Q2_correlation_%s_%s.png' %(method,config)))
plt.savefig(os.path.join(args.outdir, '%gon%g/minQ2=%g/Q2_correlation_%s_%gx%g_minQ2=%g.png' %(k,p,minq2,method,k,p,minq2)))
##########################################################################################
#function to construct Bjorken-x correlation plots #function to construct Bjorken-x correlation plots
##########################################################################################
def Xcorrelation(minq2,method): #minq2 can be 1,10,100, or 1000; method can be 'e','DA', or 'JB' def Xcorrelation(minq2,method): #minq2 can be 1,10,100, or 1000; method can be 'e','DA', or 'JB'
Xvalues_Y = method_Xvalues_dict['{}'.format(method)] #x values of the given method, that are mapped onto the y axis Xvalues_Y = method_Xvalues_dict['{}'.format(method)] #x values of the given method, that are mapped onto the y axis
...@@ -142,8 +147,7 @@ def Xcorrelation(minq2,method): #minq2 can be 1,10,100, or 1000; method can be ' ...@@ -142,8 +147,7 @@ def Xcorrelation(minq2,method): #minq2 can be 1,10,100, or 1000; method can be '
plt.xlabel('x Truth') plt.xlabel('x Truth')
plt.ylabel('$x$ {}'.format(method_dict['{}'.format(method)])) plt.ylabel('$x$ {}'.format(method_dict['{}'.format(method)]))
plt.title('{} $x$ correlation {}x{} $minQ^2=${}$GeV^2$'.format(method_dict['{}'.format(method)],k,p,minq2)) plt.title('{} $x$ correlation {}x{} $minQ^2=${}$GeV^2$'.format(method_dict['{}'.format(method)],k,p,minq2))
plt.show() plt.savefig(os.path.join(r_path, 'x_correlation_%s_%s.png' %(method,config)))
plt.savefig(os.path.join(args.outdir, '%gon%g/minQ2=%g/x_correlation_%s_%gx%g_minQ2=%g.png' %(k,p,minq2,method,k,p,minq2)))
...@@ -155,19 +159,27 @@ keys = ur.concatenate(rec_file + ':events/' + 'InclusiveKinematicsDA') ...@@ -155,19 +159,27 @@ keys = ur.concatenate(rec_file + ':events/' + 'InclusiveKinematicsDA')
DoubleAngle = [keys['InclusiveKinematicsDA.Q2'], keys['InclusiveKinematicsDA.x']] DoubleAngle = [keys['InclusiveKinematicsDA.Q2'], keys['InclusiveKinematicsDA.x']]
keys = ur.concatenate(rec_file + ':events/' + 'InclusiveKinematicsJB') keys = ur.concatenate(rec_file + ':events/' + 'InclusiveKinematicsJB')
JacquetBlondel = [keys['InclusiveKinematicsJB.Q2'], keys['InclusiveKinematicsJB.x']] JacquetBlondel = [keys['InclusiveKinematicsJB.Q2'], keys['InclusiveKinematicsJB.x']]
keys = ur.concatenate(rec_file + ':events/' + 'InclusiveKinematicsSigma')
Sigma = [keys['InclusiveKinematicsSigma.Q2'], keys['InclusiveKinematicsSigma.x']]
keys = ur.concatenate(rec_file + ':events/' + 'InclusiveKinematicseSigma')
eSigma = [keys['InclusiveKinematicseSigma.Q2'], keys['InclusiveKinematicseSigma.x']]
Q2values_T = Truth[0] Q2values_T = Truth[0]
Q2values_E = Electron[0] Q2values_E = Electron[0]
Q2values_DA = DoubleAngle[0] Q2values_DA = DoubleAngle[0]
Q2values_JB = JacquetBlondel[0] Q2values_JB = JacquetBlondel[0]
Q2values_S = Sigma[0]
Q2values_eS = eSigma[0]
Xvalues_T = Truth[1] Xvalues_T = Truth[1]
Xvalues_E = Electron[1] Xvalues_E = Electron[1]
Xvalues_DA = DoubleAngle[1] Xvalues_DA = DoubleAngle[1]
Xvalues_JB = JacquetBlondel[1] Xvalues_JB = JacquetBlondel[1]
Xvalues_S = Sigma[1]
Xvalues_eS = eSigma[1]
method_dict = {'e':'Electron','DA':'Double-Angle','JB':'Jacquet-Blondel'} method_dict = {'e':'Electron','DA':'Double-Angle','JB':'Jacquet-Blondel','S':'Sigma','eS':'eSigma'}
method_Q2values_dict = {'e':Q2values_E,'DA':Q2values_DA,'JB':Q2values_JB} method_Q2values_dict = {'e':Q2values_E,'DA':Q2values_DA,'JB':Q2values_JB,'S':Q2values_S,'eS':Q2values_eS}
method_Xvalues_dict = {'e':Xvalues_E,'DA':Xvalues_DA,'JB':Xvalues_JB} method_Xvalues_dict = {'e':Xvalues_E,'DA':Xvalues_DA,'JB':Xvalues_JB,'S':Xvalues_S,'eS':Xvalues_eS}
Q2correlation(minq2,'e') Q2correlation(minq2,'e')
Xcorrelation(minq2,'e') Xcorrelation(minq2,'e')
...@@ -175,3 +187,43 @@ Q2correlation(minq2,'DA') ...@@ -175,3 +187,43 @@ Q2correlation(minq2,'DA')
Xcorrelation(minq2,'DA') Xcorrelation(minq2,'DA')
Q2correlation(minq2,'JB') Q2correlation(minq2,'JB')
Xcorrelation(minq2,'JB') Xcorrelation(minq2,'JB')
Q2correlation(minq2,'S')
Xcorrelation(minq2,'S')
Q2correlation(minq2,'eS')
Xcorrelation(minq2,'eS')
##########################################################################################
#Kinematic coverage plot
##########################################################################################
fig = plt.figure()
gs = fig.add_gridspec(2, 3, wspace=0.09, hspace = 0.2)
(ax1, ax2, ax3), (ax4, ax5, ax6) = gs.subplots()
ax1.hist2d(ak.to_numpy(ak.flatten(Xvalues_T)), ak.to_numpy(ak.flatten(Q2values_T)), bins = [x_bins, Q2_bins])
ax2.hist2d(ak.to_numpy(ak.flatten(Xvalues_S)), ak.to_numpy(ak.flatten(Q2values_S)), bins = [x_bins, Q2_bins])
ax3.hist2d(ak.to_numpy(ak.flatten(Xvalues_DA)), ak.to_numpy(ak.flatten(Q2values_DA)), bins = [x_bins, Q2_bins])
ax4.hist2d(ak.to_numpy(ak.flatten(Xvalues_E)), ak.to_numpy(ak.flatten(Q2values_E)), bins = [x_bins, Q2_bins])
ax5.hist2d(ak.to_numpy(ak.flatten(Xvalues_eS)), ak.to_numpy(ak.flatten(Q2values_eS)), bins = [x_bins, Q2_bins])
ax6.hist2d(ak.to_numpy(ak.flatten(Xvalues_JB)), ak.to_numpy(ak.flatten(Q2values_JB)), bins = [x_bins, Q2_bins])
for axes in [ax1,ax2,ax3,ax4,ax5,ax6]:
axes.set_xscale('log')
axes.set_yscale('log')
axes.set_ylabel('$Q^2$ [$GeV^2$]')
axes.set_xlabel('$x$')
ax1.set_title('Truth')
ax2.set_title('Sigma')
ax3.set_title('Double Angle')
ax4.set_title('Electron')
ax5.set_title('eSigma')
ax6.set_title('Jacquet Blondel')
fig.set_figwidth(30)
fig.set_figheight(20)
plt.title('Kinematic Coverage {}x{} $minQ^2=${}$GeV^2$'.format(k,p,minq2))
plt.savefig(os.path.join(r_path, 'kinematic_coverage_%s.png' %(method,config)))
benchmarks/dis/dis.sh
+ 1
− 1
View file @ a93bf3b1
...@@ -122,7 +122,7 @@ if [[ "$?" -ne "0" ]] ; then ...@@ -122,7 +122,7 @@ if [[ "$?" -ne "0" ]] ; then
exit 1 exit 1
fi fi
python benchmarks/dis/analysis/kinematics_correlations.py --rec_file ${REC_FILE} --ebeam ${EBEAM} --pbeam ${PBEAM} --minq2 ${MINQ2} python benchmarks/dis/analysis/kinematics_correlations.py --rec_file ${REC_FILE} --config ${PLOT_TAG}_${DETECTOR_CONFIG} --results_path ${RESULTS_PATH}
if [[ "$?" -ne "0" ]] ; then if [[ "$?" -ne "0" ]] ; then
echo "ERROR running kinematics_correlations script" echo "ERROR running kinematics_correlations script"
exit 1 exit 1
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment