SNEWS2 / snewpy

A Python package for working with supernova neutrinos
https://snewpy.readthedocs.io
BSD 3-Clause "New" or "Revised" License
24 stars 17 forks source link

Snewpy not finding Smearing and efficiency files #270

Closed yveskini closed 9 months ago

yveskini commented 9 months ago

Hello, I am encountering an issue while attempting to run the SNOwGLoBES_usage notebook, which is located at doc/nb/SNOwGLoBES_usage.ipynb. Despite ensuring that the paths are correctly configured, I am unable to locate the smearing files and efficiency files in SNOwGLoBES. Consequently, Snewpy sets these values to 100%. To address this problem, I have made slight modifications to the _load_smearing_matrices and '_load_efficiency_vectors' methods as outlined below. As a result, I now obtain counts that align with the tutorial. I am uncertain if this is the correct approach to resolve the issue. Thanks.

def _load_efficiency_vectors(self, path): result = {} path=Path(path) for detector in self.detectors: resdet = {} ` #for file in (path/str(detector)).glob(f'effic{detector}.dat'): for file in path.glob(f'effic{detector}.dat'):` channel = file.stem[len('effic'):-len(detector)-1] logger.debug(f'Reading file ({detector},{channel}): {file}') with open(file) as f: effs = np.fromiter(f.readlines()[0].split("{")[-1].split("}")[0].split(","), float) res_det[channel]= effs result[detector]=res_det self.efficiencies = result logger.info(f'read efficiencies for detectors: {list(self.efficiencies.keys())}') logger.debug(f'efficiencies: {self.efficiencies}')

def _load_smearing_matrices(self, path): result = {} path = Path(path) for detector in self.detectors: resdet = {} `#for file in (path/str(detector)).glob(f'smear*{detector}.dat'): for file in path.glob(f'smear*{detector}.dat'):` channel = file.stem[len('smear_'):-len(detector)-1] with open(file) as f: lines = f.readlines()[1:-1] while not "{" in lines[0]: lines = lines[1:] while not "{" in lines[-1]: lines = lines[:-1] matrix = np.zeros((len(lines),len(lines))) for i,l in enumerate(lines): elements = np.fromiter(l.split("{")[-1].split("}")[0].split(","), float) matrix[i, int(elements[0]+0.1):int(elements[1]+0.1)+1] = elements[2:] res_det[channel]= matrix result[detector]=res_det self.smearings = result logger.info(f'read smearing matrices for detectors: {list(self.smearings.keys())}') logger.debug(f'smearing matrices: {self.smearings}')

Sheshuk commented 9 months ago

Hi, thank you for reporting this! Could you please provide more details? In particular: 1) What version of SNEWPY do you use? 2) What exactly are the error messages you get 3) Did you run a standard doc/nb/SNOwGLoBES_usage.ipynb, or have some modifications in it (apart from setting the SNEWPY_models_base variable to point to your models)? 4) Did you use snowglobes_data package, installed with SNEWPY (SNOwGLoBES_path = None in the notebook), or do you use your own installation of SNOwGLoBES, and tell SNEWPY to use its data via SNOwGLoBES_path?

yveskini commented 9 months ago

Thanks for your reply.

  1. I'm using SNEWPY version 1.4b2

  2. I don't really get an error, but I get a warning. I get a bunch of these warnings: /home/yves/anaconda3/envs/snewpy/lib/python3.8/site-packages/snewpy/rate_calculator.py:155: UserWarning: Smearing not found for detector=icecube, channel=ibd. Using unsmeared spectrum warn(f'Smearing not found for detector={detector}, channel={channel.name}. Using unsmeared spectrum') /home/yves/anaconda3/envs/snewpy/lib/python3.8/site-packages/snewpy/rate_calculator.py:160: UserWarning: Efficiency not found for detector=icecube, channel=ibd. Using 100% efficiency

  3. I'm running the standard doc/nb/SNOwGLoBES_usage.ipynb notebook

  4. Setting SNOwGLoBES_path = None or SNOwGLoBES_path = "/home/yves/snowglobes" gives me the same warning.

FYI, I'm a newbie and it's the first time I'm playing around with SNEWPY and SnowGlobes. Thanks

Sheshuk commented 9 months ago

I see. Can you also verify that you have the efficiency files, in your SNOwGLoBES dir:

ls /home/yves/snowglobes/icecube/

and/or in the snowglobes-data package (run this in the notebook):

import snowglobes_data
from importlib.resources import files
fpath = files(snowglobes_data)/'effic'/'icecube'
!ls $fpath

Actually I'm surprised that your patch works for you, because you try to read the efficiencies and smearing directly from the {path}, but in snowglobes they are located in {path}/{detector_name}...

yveskini commented 9 months ago

All the efficiency (smearing) files are located in "/home/yves/snowglobes/effic/" ("/home/yves/snowglobes/smear/").

ls "/home/yves/snowglobes/effic/" gives me :

effic_ibd_he_wc100kt30prct_he.dat effic_ibd_icecube.dat effic_ibd_wc100kt15prct.dat effic_ibd_wc100kt30prct.dat effic_nc_ah_nue_Ar40_ar40kt.dat effic_nc_ah_nue_Ar40_he_ar40kt_he.dat effic_nc_ah_nuebar_Ar40_ar40kt.dat effic_nc_ah_nuebar_Ar40_he_ar40kt_he.dat effic_nc_ah_numu_Ar40_ar40kt.dat effic_nc_ah_numu_Ar40_he_ar40kt_he.dat effic_nc_ah_numubar_Ar40_ar40kt.dat effic_nc_ah_numubar_Ar40_he_ar40kt_he.dat effic_nc_ah_nutau_Ar40_ar40kt.dat effic_nc_ah_nutau_Ar40_he_ar40kt_he.dat effic_nc_ah_nutaubar_Ar40_ar40kt.dat effic_nc_ah_nutaubar_Ar40_he_ar40kt_he.dat effic_nc_nue_Ar40_ar40kt.dat effic_nc_nue_Ar40_he_ar40kt_he.dat effic_nc_nuebar_Ar40_ar40kt.dat effic_nc_nuebar_Ar40_he_ar40kt_he.dat effic_nc_nuebar_O16_icecube.dat effic_nc_nuebar_O16_wc100kt15prct.dat effic_nc_nuebar_O16_wc100kt30prct.dat effic_nc_nuebar_O16_wc100kt30prct_he.dat effic_nc_nuebar_Pb208_1n_halo1.dat effic_nc_nuebar_Pb208_1n_halo2.dat effic_nc_nuebar_Pb208_2n_halo1.dat effic_nc_nuebar_Pb208_2n_halo2.dat effic_nc_nue_O16_icecube.dat effic_nc_nue_O16_wc100kt15prct.dat effic_nc_nue_O16_wc100kt30prct.dat effic_nc_nue_O16_wc100kt30prct_he.dat effic_nc_nue_Pb208_1n_halo1.dat effic_nc_nue_Pb208_1n_halo2.dat effic_nc_nue_Pb208_2n_halo1.dat effic_nc_nue_Pb208_2n_halo2.dat effic_nc_numu_Ar40_ar40kt.dat effic_nc_numu_Ar40_he_ar40kt_he.dat effic_nc_numubar_Ar40_ar40kt.dat effic_nc_numubar_Ar40_he_ar40kt_he.dat effic_nc_numubar_O16_icecube.dat effic_nc_numubar_O16_wc100kt15prct.dat effic_nc_numubar_O16_wc100kt30prct.dat effic_nc_numubar_O16_wc100kt30prct_he.dat effic_nc_numubar_Pb208_1n_halo1.dat effic_nc_numubar_Pb208_1n_halo2.dat effic_nc_numubar_Pb208_2n_halo1.dat effic_nc_numubar_Pb208_2n_halo2.dat effic_nc_numu_O16_icecube.dat effic_nc_numu_O16_wc100kt15prct.dat effic_nc_numu_O16_wc100kt30prct.dat effic_nc_numu_O16_wc100kt30prct_he.dat effic_nc_numu_Pb208_1n_halo1.dat effic_nc_numu_Pb208_1n_halo2.dat effic_nc_numu_Pb208_2n_halo1.dat effic_nc_numu_Pb208_2n_halo2.dat effic_nc_nutau_Ar40_ar40kt.dat effic_nc_nutau_Ar40_he_ar40kt_he.dat effic_nc_nutaubar_Ar40_ar40kt.dat effic_nc_nutaubar_Ar40_he_ar40kt_he.dat effic_nc_nutaubar_O16_icecube.dat effic_nc_nutaubar_O16_wc100kt15prct.dat effic_nc_nutaubar_O16_wc100kt30prct.dat effic_nc_nutaubar_O16_wc100kt30prct_he.dat effic_nc_nutaubar_Pb208_1n_halo1.dat effic_nc_nutaubar_Pb208_1n_halo2.dat effic_nc_nutaubar_Pb208_2n_halo1.dat effic_nc_nutaubar_Pb208_2n_halo2.dat effic_nc_nutau_O16_icecube.dat effic_nc_nutau_O16_wc100kt15prct.dat effic_nc_nutau_O16_wc100kt30prct.dat effic_nc_nutau_O16_wc100kt30prct_he.dat effic_nc_nutau_Pb208_1n_halo1.dat effic_nc_nutau_Pb208_1n_halo2.dat effic_nc_nutau_Pb208_2n_halo1.dat effic_nc_nutau_Pb208_2n_halo2.dat effic_nue_Ar40_ar40kt.dat effic_nue_Ar40_he_ar40kt_he.dat effic_nue_Ar40_klmv_ar40kt.dat effic_nue_Ar40_marley1_ar40kt.dat effic_nue_Ar40_marley2_ar40kt.dat effic_nuebar_Ar40_ar40kt.dat effic_nuebar_Ar40_he_ar40kt_he.dat effic_nuebar_Ar40_klmv_ar40kt.dat effic_nuebar_e_ar40kt.dat effic_nuebar_e_halo1.dat effic_nuebar_e_halo2.dat effic_nuebar_e_he_ar40kt_he.dat effic_nuebar_e_he_wc100kt30prct_he.dat effic_nuebar_e_icecube.dat effic_nuebar_e_wc100kt15prct.dat effic_nuebar_e_wc100kt30prct.dat effic_nuebar_O16_icecube.dat effic_nuebar_O16_wc100kt15prct.dat effic_nuebar_O16_wc100kt30prct.dat effic_nuebar_O16_wc100kt30prct_he.dat effic_nue_e_ar40kt.dat effic_nue_e_halo1.dat effic_nue_e_halo2.dat effic_nue_e_he_ar40kt_he.dat effic_nue_e_he_wc100kt30prct_he.dat effic_nue_e_icecube.dat effic_nue_e_wc100kt15prct.dat effic_nue_e_wc100kt30prct.dat effic_nue_O16_icecube.dat effic_nue_O16_wc100kt15prct.dat effic_nue_O16_wc100kt30prct.dat effic_nue_O16_wc100kt30prct_he.dat effic_nue_Pb208_1n_halo1.dat effic_nue_Pb208_1n_halo2.dat effic_nue_Pb208_2n_halo1.dat effic_nue_Pb208_2n_halo2.dat effic_numubar_e_ar40kt.dat effic_numubar_e_halo1.dat effic_numubar_e_halo2.dat effic_numubar_e_he_ar40kt_he.dat effic_numubar_e_he_wc100kt30prct_he.dat effic_numubar_e_icecube.dat effic_numubar_e_wc100kt15prct.dat effic_numubar_e_wc100kt30prct.dat effic_numubar_H1_wc100kt30prct_he.dat effic_numubar_O16_wc100kt30prct_he.dat effic_numu_e_ar40kt.dat effic_numu_e_halo1.dat effic_numu_e_halo2.dat effic_numu_e_he_ar40kt_he.dat effic_numu_e_he_wc100kt30prct_he.dat effic_numu_e_icecube.dat effic_numu_e_wc100kt15prct.dat effic_numu_e_wc100kt30prct.dat effic_numu_O16_wc100kt30prct_he.dat effic_nutaubar_e_ar40kt.dat effic_nutaubar_e_halo1.dat effic_nutaubar_e_halo2.dat effic_nutaubar_e_he_ar40kt_he.dat effic_nutaubar_e_he_wc100kt30prct_he.dat effic_nutaubar_e_icecube.dat effic_nutaubar_e_wc100kt15prct.dat effic_nutaubar_e_wc100kt30prct.dat effic_nutau_e_ar40kt.dat effic_nutau_e_halo1.dat effic_nutau_e_halo2.dat effic_nutau_e_he_ar40kt_he.dat effic_nutau_e_he_wc100kt30prct_he.dat effic_nutau_e_icecube.dat effic_nutau_e_wc100kt15prct.dat effic_nutau_e_wc100kt30prct.dat

I also did check that I'm using the correct version of Snowglobes ( V1.2)

import snowglobes_data from importlib_resources import files fpath = files(snowglobes_data)/'effic'/'icecube' !ls $fpath gives me this:

effic_ibd_icecube.dat effic_nuebar_O16_icecube.dat effic_nc_nuebar_O16_icecube.dat effic_nue_e_icecube.dat effic_nc_nue_O16_icecube.dat effic_nue_O16_icecube.dat effic_nc_numubar_O16_icecube.dat effic_numubar_e_icecube.dat effic_nc_numu_O16_icecube.dat effic_numu_e_icecube.dat effic_nc_nutaubar_O16_icecube.dat effic_nutaubar_e_icecube.dat effic_nc_nutau_O16_icecube.dat effic_nutau_e_icecube.dat effic_nuebar_e_icecube.dat

Without my patch (using Snewpy as it is), doing this instead (below) fixes my problem: import snowglobes_data SNOwGLoBES_path = files(snowglobes_data)

Thanks a lot

Sheshuk commented 9 months ago

Oh I see the problem: in SNOwGLoBES v1.3 the directory structure was changed compared to 1.2, so this leads to problems like this. I see that our current documentation for "stable" release still recommends v1.2, which is not right. We'll need to update this, thanks again for reporting the problem!

Still, just leaving SNOwGLoBES_path = None in the example notebook should lead to SNEWPY using snowglobes_data package internally. Unless you have $SNOWGLOBES environment var set to your local snowglobes v1.2 installation. So I think the solution should be just to unset SNOWGLOBES before running the notebook.

yveskini commented 9 months ago

Nice, thank you!

If I may also ask, is the 'Garching' model from Huedepohl et al. (http://link.aps.org/doi/10.1103/PhysRevLett.104.251101), with the fluxes provided by default in 'Snoglobes,' also a default model in SNEWPY? If not, is it something that one can provide to SNEWPY without modifying or hacking it to match? This model should be close enough to the Bollig_2016 (defined ccsn.py).

JostMigenda commented 9 months ago

@Sheshuk Unfortunately, the current stable release of snewpy (v1.3) still expects SNOwGLoBES v1.2, so the documentation is correct there. But issues with incompatible snewpy/SNOwGLoBES versions pop up regularly, so that should be a reminder to us to get snewpy v1.4 out as soon as possible. I’ll look into it Monday.

@yveskini The Hüdepohl model is not officially supported by snewpy right now, no. (We try to keep up with modern simulations and are less focussed on adding older ones; though you’re right that it might be worth making an exception for the SNOwGLoBES standard model.) For now, if the file format is identical to the one used by Bollig_2016, you can use the following little hack:

from snewpy.models.loaders import Bollig_2016
from os.path import abspath
model = Bollig_2016(abspath('path/to/Huedepohl_2010/s15.0c_3D_fastrot_dir1')) # drop the `_LS220_nue` ending of the filename

… and then the example code in models/Bollig_2016/Bollig_2016.ipynb should work.

To get SNOwGLoBES_usage.ipynb to run will require a few tweaks to snewpy itself; off the top of my head I think adding

class Huedepohl_2010(Bollig_2016):
    pass

in ccsn.py should be sufficient?

If there are any changes in the input file format, things will be more complicated. In that case I’d double check whether you really need to use that particular model instead of the more recent ones already available in snewpy; but if the answer is yes, then we’ll help you with that.

Sheshuk commented 9 months ago

Unfortunately, the current stable release of snewpy (v1.3) still expects SNOwGLoBES v1.2, so the documentation is correct there.

You're right. Probably for it to be less misleading, we should add the explicit SNEWPY version in the documentation in header or some other visible place. Otherwise a user installs some release and opens "stable" documentation version, and can't crosscheck if it is relevant or not.

yveskini commented 9 months ago

@JostMigenda Thanks a lot for the tips! My naive initial approach was to create input files for each flavor in a format similar to the Bollig model, using the fluxes provided in Snoglobes. However, upon closer inspection of the Snoglobes flux data, it became apparent that this task isn't as straightforward as I initially thought. The timestamps for luminosity and average energies don't necessarily align.

Nevertheless, I appreciate your assistance. I'll keep an eye out for any further suggestions or guidance.