nu-radio / NuRadioReco

reconstruction framework for radio detectors of high-energy neutrinos
GNU General Public License v3.0
5 stars 3 forks source link

Galactic noise simulation #269

Closed christophwelling closed 3 years ago

christophwelling commented 3 years ago

Adds a module to simulate noise from radio sources in the sky. It uses the PyGDSM package to get a radio map of the sky, evaluates it at a set of points to get the radio emission from different directions and folds it with the antenna response.

christophwelling commented 3 years ago

Good idea. The reason I did it this way was that pygdsm can provide the map in horizon coordinates, so this was more convenient. But that seemed a bit bugged, so I implemented it by hand. But you are right, now that I do it this way there is no reason to stick to that.

christophwelling commented 3 years ago

Okay, I now downsize the skymap to have a manageable number of pixels and then simulate the noise emission for each one of them (unless they are below the horizon).

anelles commented 3 years ago

In LOFAR we never did the calculations in the time domain, integrated power only. But I think we need the phases to be as random as possible.

On Fri, Nov 6, 2020, 10:13 Christian Glaser notifications@github.com wrote:

@cg-laser commented on this pull request.

In NuRadioReco/modules/channelGalacticNoiseAdder.py https://github.com/nu-radio/NuRadioReco/pull/269#discussion_r518617123:

  • temperature_interpolator = scipy.interpolate.interp1d(self.__interpolaiton_frequencies, np.log10(noise_temperatures[:, i_zenith, i_azimuth]), kind='quadratic')
  • noise_temperature = np.power(10, temperature_interpolator(freqs[passband_filter]))
  • S = (2. scipy.constants.Boltzmann (freqs[passband_filter] / units.Hz)2 / scipy.constants.c2 (noise_temperature / units.kelvin) solid_angle) * (units.watt / units.m**2 / units.Hz)
  • S[np.isnan(S)] = 0
  • S_per_bin = S * d_f
  • flux_sum += S_per_bin
  • E = np.sqrt(S_per_bin / (units.watt / units.m*2) / (scipy.constants.c scipy.constants.epsilon_0)) / (d_f) * units.V / units.m
  • if self.__debug:
  • ax_1.scatter(self.__interpolaiton_frequencies / units.MHz, noise_temperatures[:, i_zenith, i_azimuth] / units.kelvin, c='k', alpha=.01)
  • ax_1.plot(freqs[passband_filter] / units.MHz, noise_temperature, c='k', alpha=.02)
  • ax_2.plot(freqs[passband_filter] / units.MHz, S_per_bin / d_f / (units.watt / units.m**2 / units.MHz), c='k', alpha=.02)
  • ax_3.plot(freqs[passband_filter] / units.MHz, E / (units.V / units.m), c='k', alpha=.02)
  • noise_spectrum = np.zeros((3, freqs.shape[0]), dtype=np.complex)
  • phases = np.random.uniform(0, 2. * np.pi, len(S))
  • polarizations = np.random.uniform(0, 2. * np.pi, len(S))

Maybe @anelles https://github.com/anelles remembers how it was done in the Lofar analysis? Maybe instead the phi and theta component get different random phases?

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/nu-radio/NuRadioReco/pull/269#discussion_r518617123, or unsubscribe https://github.com/notifications/unsubscribe-auth/AC663BOGDRLGCB6M55SXCETSOO45BANCNFSM4SRXZZTA .