Example: Full reconstruction cycle

This example shows a full reconstruction cycle for simulated CoREAS events. The full cycle consists of event generation, event selection and event reconstruction. The full python code of this example is available (HERE). For running the reconstruction, just execute

python FullReconstruction.py 32 test_data.hdf5

CoREAS reader and event writer

Several CoREAS reader modules are provided, reading the hdf5 simulation files of CoREAS events. In this example, the readCoREAS module is used, which will generate a number of core positions randomly distributed and it determines the closest simulated station and the corresponding event object. In this way, a realistic distribution of cosmic-ray events is obtained and it allows to repeatedly use the same CoREAS simulation.

readCoREAS.begin(input_files, station_id, n_cores=10, max_distance=None)
for iE, evt in enumerate(readCoREAS.run(det)):

For each event in the CoREAS files the reconstruction cycle is run. All the data of events that pass the selection are written in NuRadioReco format by the event writer at the end of the cycle, by:

eventWriter.run(evt)

Event generation

From the CoREAS files, the events are selected which have a signal in the frequency band of 100 - 500 MHz, meaning an amplitude of at least 8 times the standard deviation of the noise. The voltage traces produced in each channel by the electric-field are calculated and the resulting traces are written into the channel storage. To fully simulate data-like events, the phase and gain of the ARIANNA hardware are taken into account and noise is added to the voltage traces. A filter is applied to the channels and the simulated electric-field is downsampled in order to reduce file size.

The corresponding code is given by:

if simulationSelector.run(evt, station.get_sim_station(), det):
   efieldToVoltageConverter.run(evt, station, det)
   hardwareResponseIncorporator.run(evt, station, det, sim_to_data=True)
   channelGenericNoiseAdder.run(evt, station, det, type = "perfect_white", amplitude = 20* units.mV)
   electricFieldResampler.run(evt, station.get_sim_station(), det, sampling_rate=1 * units.GHz)

Event selection

For the event selection a trigger is simulated, i.e., if a certain event is recorded by the detector. Only the events marked as triggered are selected. The events are filtered in a frequency range of 50-800 MHz.

The corresponding code is given below:

triggerSimulator.run(evt, station,det, number_concidences = 2, threshold = 100 *units.mV)
if station.get_trigger('default_simple_threshold').has_triggered():
   channelBandPassFilter.run(evt, station, det, passband=[80 * units.MHz, 500 * units.MHz], filter_type='butter', order = 10)

Event reconstruction

An additional filter is applied to the channel traces, to get rid of the amplitudes in the small frequencies. In order to reconstruct the fully simulated and selected events, there is corrected for the ARIANNA hardware. The direction of the event is fitted and used for the reconstruction of the electric-field. The electric-field is reconstructed using the forward folding technique, which fits an analytic model of the electric-field pulse to voltage traces. The voltage- and electric-field traces are downsampled to reduce the file size.

channelBandPassFilter.run(evt, station, det, passband=[60 * units.MHz,600 * units.MHz], filter_type='rectangular')
hardwareResponseIncorporator.run(evt, station, det)
templateDirectionFitter.run(evt, station, det, cosmic_ray=True, channels_to_use=used_channels_fit)
correlationDirectionFitter.run(evt, station, det, n_index=1., channel_pairs=channel_pairs)
voltageToAnalyticEfieldConverter.run(evt, station, det, use_channels=used_channels_efield, bandpass=[80*units.MHz, 500*units.MHz], useMCdirection=False)
channelResampler.run(evt, station, det, sampling_rate=1 * units.GHz)
electricFieldResampler.run(evt, station, det, sampling_rate=1 * units.GHz)