ADD overview_sums.py

AHBGOx.py

@@ -20,6 +20,7 @@ parser.add_argument("-two", action='store_true', help = "if chosen, functions fu
 parser.add_argument("-zero", action='store_true', help = "if chosen, data will be normalized to 0°C, chose only for sum")  # on/off flag
 parser.add_argument("-export", type=str, nargs='+', help = "filpath to store output if other then currend directory/histogramms/filename")
 parser.add_argument("-nameadd", type=str, nargs='+', help = "custom nameadd to add to filename")
+parser.add_argument("-six", action='store_true', help = "use channels for 'AHBGOS'")
 #parser.add_argument("-tdhist",  action='store_true', help = "makes 2d histogramm of B1 against B2 if AHBGOC")
 
 args = parser.parse_args()
@@ -47,6 +48,8 @@ if 'AHBGOD' in file:
     x = 'D'
 elif 'AHBGOC' in file:
     x = 'C'
+    if args.six:
+        x='S'
 elif 'AHBGOS' in file:
     x = 'S'
 

overview.py

@@ -43,6 +43,7 @@ layoutchange=False
 if args.filenumber == 19 or args.filenumber ==20:
     layoutchange = True
     fpos = [2,5,7,12,13,16]
+    fpos = [16]
 filename = "histograms/" + file + "_B1B2sum_vertical"    
 
 def main(filename, fpos):
@@ -68,6 +69,8 @@ def main(filename, fpos):
         #num  [0,1,2,3,4,5,6, 7, 8,9,10,11,12,13,14,15,16]
 
         #channels=['B1','B2','sum']
+        
+#Define startparameters for fit
     a = [165,165,220]
     b = [290,290,400]
     startparams = [[2400,65,7,160,60],[2400,65,7,160,60],[1300,115,15,80,130]]
@@ -75,7 +78,7 @@ def main(filename, fpos):
     if layoutchange:
         a = 6*[140]+[220]
         b = 6*[175]+[400]
-        startparams = 6*[[2400,25,3,160,60]] + [[1300,115,15,80,130]]
+        startparams = 6*[[2400,26,3,160,60]] + [[1300,115,15,80,130]]
         
     resV = 0.838214/2
 
@@ -119,7 +122,7 @@ def main(filename, fpos):
                     par, cov = curve_fit(lan, data[energy][a[c]:b[c]], data[channels[c]][a[c]:b[c]], startparams[c], sigma = sig[a[c]:b[c]], absolute_sigma=True)
                     print(par)
                     perr = np.sqrt(np.diag(cov))
-                    print(perr)
+                    #print(perr)
                     axs[i].step(np.multiply(data[energy],1.0/f), np.multiply(data[channels[c]], 3600/resV/Itime*f), label=channels[c]+' ; $e$ = ('+ str(round(par[1]/f,2)) + '$\pm$' +str(round(perr[1]/f+0.005,2))+') '+energystr)
                     axs[i].plot(np.multiply(data[energy][a[c]:b[c]],1.0/f), np.multiply(lan(data[energy], *par), 3600/resV/Itime*f)[a[c]:b[c]], 'k--') 
               

overview_sums.py

@@ -0,0 +1,147 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Created on Thu May  2 13:21:02 2024
+
+@author: ava
+"""
+import matplotlib.pyplot as plt
+import argparse
+import sys
+import numpy as np
+import pandas as pd
+import os.path
+from scipy.optimize import curve_fit
+from fit_functions import lan
+
+def main():
+
+    parser = argparse.ArgumentParser(description="This script plots overview for position dependencies in BigBGO")
+    parser.add_argument("-v", "--verbose", choices = ['0','1','2'], help = "when zero only emit serious error messages; when one also emit warnings, when two emit informational messages")
+    parser.add_argument("filenumbers", type = int,nargs='+', help="number of files of wich data shall be plottet (14, 16, 18,19,20)")
+    parser.add_argument("-f", "--fit", action='store_true', help = "landau distributions will be fittet to data. Has to be adjusted for 19!!!")  # on/off flag
+    parser.add_argument("-p", "--photoelectrons", action='store_true', help = "converts data from keV{Si} to amount of photoelectrons")
+    parser.add_argument("-e", "--export", type = str, nargs='+', help = "first: filepath were figure shall be exported to (if f, filename will be 'overview-[number]') \n second:format (png,pdf,...)")
+
+    args = parser.parse_args()
+
+    files=[]
+    numberstr=''
+    for num in args.filenumbers:
+        numberstr = numberstr + '-'+str(num)
+        if num == 14:
+            files.append("2023-12-21-AHBGOC-14-langzeit")
+        if num == 16:
+            files.append("2024-01-19-AHBGOC-rau-16_2weeks")
+        if num == 18:
+            files.append("2024-03-27-AHBGOC-six-18")
+        if num == 19:
+            files.append("2024-04-08-AHBGOC-six_30V-19")
+        if num == 20:
+            files.append("2024-04-11-AHBGOC-six-20")
+    
+    if len(files)<1:
+        print("No valid filenamenumber, choose from '14', '16', '18' , '19' ,'20'.")
+        sys.exit()
+        
+    fpos = [2,5,7,9,11,12,13,16]
+
+    energy = 'keV'
+    if args.photoelectrons:
+        f = 0.00362
+        energystr = '$e_{ph}$'
+    else: 
+        f = 1.0
+        energystr = energy + '$_{Si}$'
+        
+        
+    plt.figure(figsize=(15,10))    
+    process_data(files,f,energy,fpos,args.fit,args.filenumbers)
+    if len(args.filenumbers)>1:
+        tit = "Overview sums of files " + numberstr
+    else:
+        tit = "Overview sums of '"+files[0]+"'"
+    layout(tit,energystr,f)
+    
+    if not args.export:
+        plt.show()
+    else:
+        if args.export[0] == 'f':
+            name = 'plots/overview_sums'+numberstr
+        else:
+            name = args.export[0]
+        if len(args.export)>1:
+                end = args.export[1]
+        else:
+                end=''
+        plt.savefig(name+'.'+end)
+        
+        
+ ################################################################       
+        
+        
+        
+    
+def read_Itime(file):
+  
+    f = file + '.Itime'
+    Itime = None 
+    while len(f) > 0:
+        if os.path.isfile(f+'.Itime'):
+            Itime = int(pd.read_csv(f+'.Itime', sep='\s+', skiprows=0, names = ['time'])['time'][0])
+            f=''
+        else:
+                f=f[:-1]
+
+    return Itime
+
+def process_data(files,f,energy,fpos,fit,fnum):
+    resV = 0.838214/2
+    colors=['b','r','g','m','c']
+    for i in range(len(files)):
+        print(files[i])
+        filename = "histograms/" + files[i] + "_B1B2sum_vertical" 
+        Itime=read_Itime(filename)  
+        for p in fpos:
+            channels=[]
+            ch = pd.read_csv(filename+"%02d" % p+".hist", sep = ' ', nrows=0)
+            for c in ch:
+                channels.append(str(c))
+            channels = channels[1:]
+        
+            data = pd.read_csv(filename+"%02d" % p+".hist", sep='\s+')
+    
+            for c in range(len(channels)):
+                if channels[c]=='sum':
+                    if not fit:
+                        plt.step(np.multiply(data[energy],1.0/f),np.multiply(data[channels[c]], 3600/resV/Itime*f), color=colors[i],label=channels[c]+'-P'+str(p)+'-'+str(fnum[i]))
+                    else:
+                        sig = np.power(data[channels[c]]+1,0.5)
+                        par, cov = curve_fit(lan, data[energy][a[c]:b[c]], data[channels[c]][a[c]:b[c]], startparams[c], sigma = sig[a[c]:b[c]], absolute_sigma=True)
+                        print(par)
+                        perr = np.sqrt(np.diag(cov))
+                        #print(perr)
+                        plt.step(np.multiply(data[energy],1.0/f), np.multiply(data[channels[c]], 3600/resV/Itime*f), label=channels[c]+' ; $e$ = ('+ str(round(par[1]/f,2)) + '$\pm$' +str(round(perr[1]/f+0.005,2))+') '+energystr)
+                        plt.plot(np.multiply(data[energy][a[c]:b[c]],1.0/f), np.multiply(lan(data[energy], *par), 3600/resV/Itime*f)[a[c]:b[c]], 'k--') 
+        
+    
+def layout(tit,energystr,f):
+
+    fs = 15
+    lim = [30/f,280/f, 0.003*f,1.8*f]
+
+    ylab = 'counts [1/h/' +energystr+']'
+    xlab = 'A [' + energystr+']'
+    plt.legend(fontsize = fs)
+    plt.xticks(fontsize=fs)
+    plt.yticks(fontsize=fs)
+    plt.title(tit, fontsize = fs+10, pad = 15)
+
+    plt.yscale('log')
+    plt.ylabel(ylab, fontsize = fs+5)
+    plt.xlabel(xlab, fontsize=fs+5)
+    plt.ylim(lim[2],lim[3])
+    plt.xlim(lim[0],lim[1])
+    
+    
+main()