diff --git a/python/convertRoot2Fits.py b/python/convertRoot2Fits.py index 090dbff732a1e834d771e5f546182b30de394f29..646de19cbbb01814c060fef08bd3401b0f6fb659 100644 --- a/python/convertRoot2Fits.py +++ b/python/convertRoot2Fits.py @@ -102,4 +102,4 @@ def generate_fits(energy, theta_deg, phi_deg): except Exception as e: print("Errore durante la creazione del file .fits") -generate_fits(66, 0, 0) \ No newline at end of file +generate_fits(66, 40, 40) \ No newline at end of file diff --git a/python/data.py b/python/data.py index 23148c0496e900e1ab2a3509065dfab5028e9035..f5c4b3b3c1f1fd71ce4807d6af0ca5ca8bd599b8 100644 --- a/python/data.py +++ b/python/data.py @@ -6,7 +6,7 @@ import astropy.io.fits as pyfits from noisy import add_noise_X, add_noise_S -def write_fits_dep(energy, en_primary, theta, phi, type, en_dep_noisy, en_dep, name_id, pos_id, x_pixel, y_pixel): # ! è qui che decido quali informazioni sulla simulazione voglio nei file fits_dep +def write_fits_dep(event_id, energy, en_primary, theta, phi, type, en_dep_noisy, en_dep, name_id, pos_id, x_pixel, y_pixel): # ! è qui che decido quali informazioni sulla simulazione voglio nei file fits_dep x1_pixel = (pos_id // 89) +1 y1_pixel = (pos_id % 89) +1 @@ -19,13 +19,14 @@ def write_fits_dep(energy, en_primary, theta, phi, type, en_dep_noisy, en_dep, n # Extension tbhdu = pyfits.BinTableHDU.from_columns([ - #pyfits.Column(name='En_dep_Noise', format='1D', unit='keV', array = en_dep_noisy), + pyfits.Column(name='EventID', format='1J', array = event_id), pyfits.Column(name='En_Primary', format='1D', unit='keV', array = en_primary), + pyfits.Column(name='En_dep', format='1D', unit='keV', array = en_dep), + pyfits.Column(name='En_dep_Noise', format='1D', unit='keV', array = en_dep_noisy), #pyfits.Column(name=name_id, format='1I', array = pos_id), - pyfits.Column(name='X_Y', format='2J', array = X_Y), + pyfits.Column(name='X_Y', format='2I', array = X_Y), pyfits.Column(name='X_Pixel_cm', format='1D', unit='cm', array = x_pixel), - pyfits.Column(name='Y_Pixel_cm', format='1D', unit='cm', array = y_pixel) - #pyfits.Column(name='En_dep', format='1D', unit='keV', array = en_dep) + pyfits.Column(name='Y_Pixel_cm', format='1D', unit='cm', array = y_pixel) ]) tbhdu.header['EXTNAME'] = 'EVENTS' @@ -66,35 +67,64 @@ def get_data(files): with pyfits.open(files) as hdul: data = hdul[1].data tot_events = len(data) + pixel_id_col = data['Pixel_ID'] pixel_id = pixel_id_col[pixel_id_col != -1000] + scint_id_col = data['Scint_ID'] scint_id = scint_id_col[scint_id_col != -1000] + + # mask_X = (data['Scint_ID'] == -1000) # condizione //TOdO questo è un altro modo di gestire i dati in modo piu compatto, non è concluso + # eventi_X = data[mask_X] + # mask_S = (data['Scint_ID'] != -1000) # condizione + # eventi_S = data[mask_S] + # columns_for_X = ['EventID','En_Dep','Pixel_ID','X_Pixel','Y_Pixel','En_Primary'] # qui scelgo le info + # columns_for_S = ['EventID','En_Dep','Scint_ID','X_Pixel','Y_Pixel','En_Primary'] + # eventi_X_reduced = eventi_X[columns_for_X] + # eventi_S_reduced = eventi_S[columns_for_S] + + event_id = data['EventID'] + X_event_id = event_id[scint_id_col == -1000] + S_event_id = event_id[scint_id_col != -1000] + energy = data['En_Primary'][0] + en_primary = data['En_Primary'] + X_en_primary = en_primary[scint_id_col == -1000] + S_en_primary = en_primary[scint_id_col != -1000] + en_dep = data['En_dep'] - x_pixel = data['X_Pixel'] - y_pixel = data['Y_Pixel'] X_en_dep = en_dep[scint_id_col == -1000] S_en_dep = en_dep[scint_id_col != -1000] - eventi_X = (scint_id_col == -1000).sum() - eventi_S = tot_events - eventi_X - ratio_X = eventi_X/tot_events - ratio_S = eventi_S/tot_events + + x_pixel = data['X_Pixel'] + X_x_pixel = x_pixel[scint_id_col == -1000] + S_x_pixel = x_pixel[scint_id_col != -1000] + + y_pixel = data['Y_Pixel'] + X_y_pixel = y_pixel[scint_id_col == -1000] + S_y_pixel = y_pixel[scint_id_col != -1000] + + num_eventi_X = (scint_id_col == -1000).sum() + num_eventi_S = (scint_id_col != -1000).sum() + ratio_X = num_eventi_X/tot_events + ratio_S = num_eventi_S/tot_events + if num_eventi_X + num_eventi_S != tot_events: + raise ValueError("Errore: Eventi X + Eventi S DIVERSO da Eventi Totali!") X_en_dep_noisy = [add_noise_X(X_dep) for X_dep in X_en_dep] # list comprehension di Python, una struttura compatta che consente di applicare una funzione a ciascun elemento di un array o di una lista, creando direttamente un nuovo array con i risultati. S_en_dep_noisy = [add_noise_S(S_dep) for S_dep in S_en_dep] - return energy, en_primary, ratio_X, ratio_S, theta, phi, X_en_dep_noisy, S_en_dep_noisy, X_en_dep, S_en_dep, pixel_id, scint_id, x_pixel, y_pixel + return X_event_id, S_event_id, energy, X_en_primary, S_en_primary, ratio_X, ratio_S, X_x_pixel, S_x_pixel, X_y_pixel, S_y_pixel, theta, phi, X_en_dep_noisy, S_en_dep_noisy, X_en_dep, S_en_dep, pixel_id, scint_id # scrive sia il file di testo che i file fits def write_results(file_dat, results): with open(file_dat, "w") as f: - for energy, en_primary, ratio_X, ratio_S, theta, phi, X_en_dep_noisy, S_en_dep_noisy, X_en_dep, S_en_dep, pixel_id, scint_id, x_pixel, y_pixel in results: + for X_event_id, S_event_id, energy, X_en_primary, S_en_primary, ratio_X, ratio_S, X_x_pixel, S_x_pixel, X_y_pixel, S_y_pixel, theta, phi, X_en_dep_noisy, S_en_dep_noisy, X_en_dep, S_en_dep, pixel_id, scint_id in results: - X_fits = write_fits_dep(energy, en_primary, theta, phi, 'X', X_en_dep_noisy, X_en_dep, 'Pixel_ID', pixel_id, x_pixel, y_pixel) - S_fits = write_fits_dep(energy, en_primary, theta, phi, 'S', S_en_dep_noisy, S_en_dep, 'Scint_ID', scint_id, x_pixel, y_pixel) + X_fits = write_fits_dep(X_event_id, energy, X_en_primary, theta, phi, 'X', X_en_dep_noisy, X_en_dep, 'Pixel_ID', pixel_id, X_x_pixel, X_y_pixel) + S_fits = write_fits_dep(S_event_id, energy, S_en_primary, theta, phi, 'S', S_en_dep_noisy, S_en_dep, 'Scint_ID', scint_id, S_x_pixel, S_y_pixel) #f.write(f"{energy:.1f}\t{ratio_X:.3f}\t{ratio_S:.3f}\t{theta}\t{phi}\t{X_fits}\t{S_fits}\n")