TauRunner doesn't seem to be able to calculate the column depth for certain angles.

[williamluszczak]

Certain choices of angle will result in TauRunner throwing an error when trying to calculat column depth. E.g:

    Earth    = construct_earth(layers=[(3., 1.0)])
    thetas_nadir = np.array([66.44148063158723])
    tracks = make_tracks(thetas_nadir)
    mytrack = tracks[thetas_nadir[0]]

    print("column depth", mytrack.total_column_depth(Earth))

Results in an error:

  File "generate_evts.py", line 116, in <module>
    print("column depth", mytrack.total_column_depth(Earth))
  File "/home/wluszczak/taurunner/taurunner/track/track.py", line 53, in total_column_depth
    set_spline(self, body)
  File "/home/wluszczak/taurunner/taurunner/track/utils/spline_column_depth.py", line 159, in set_spline
    X2x = construct_X2x(x2X)
  File "/home/wluszczak/taurunner/taurunner/track/utils/spline_column_depth.py", line 100, in construct_X2x
    raise ValueError('Constant X2x is only allow for tracks with no column depth')
ValueError: Constant X2x is only allow for tracks with no column depth```

[jlazar17]

TauRunner requires input in radians when running as an imported module. We have added a warning if the nadir angle is outside of [0, π].

[williamluszczak]

I still seem to be getting some problems with values in that range:

thetas_nadir = np.array([2.323])
tracks = make_tracks(thetas_nadir)
mytrack = tracks[thetas_nadir[0]]

print("column depth", mytrack.total_column_depth(Earth))

Gives the same column depth error as before, but weirdly 2.32 radians works just fine

[tbertolez]

Sorry for reopening this, but I seem to have the same problem. Using a code very similar to the ones in the examples (but using an isotropic flux over the whole sky)

nevents  = 50 # number of events to simulate
eini     = 1e19 # initial energy in eV
th_min = 0
th_max = 180 # we want an isotropic flux over the whole 4pi solid angle
pid      = 16 # PDG MC Encoding particle ID (nutau)
xs_model = "CSMS" # neutrino cross section model 

Earth    = construct_earth(layers=[(4., 1.0)]) # Make Earth object with 4km water layer
xs       = CrossSections(xs_model)
rand  = np.random.RandomState(seed=7)

energies = make_initial_e(nevents, eini) # Return array of initial energies in eV
thetas   = make_initial_thetas(nevents, (th_min, th_max), rand = rand)

tau_prop = make_propagator(pid, Earth)

output = run_MC(energies, 
                thetas, 
                Earth, 
                xs,
                tau_prop, 
                rand,
         )

I get the following error:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-2-4e57f02d4ef1> in <module>
     22                 xs,
     23                 tau_prop,
---> 24                 rand,
     25          )

[~/.local/lib/python3.7/site-packages/taurunner/main.py](https://vscode-remote+ssh-002dremote-002blogin-002erc-002efas-002eharvard-002eedu.vscode-resource.vscode-cdn.net/n/home11/tbertolez/IceCubeZoomEeV/~/.local/lib/python3.7/site-packages/taurunner/main.py) in run_MC(eini, thetas, body, xs, propagator, rand, no_secondaries, flavor, no_losses, condition, depth, track_type)
    210         )
    211 
--> 212         out      = Propagate(particle, my_track, body, condition=condition)
    213 
    214         if (out.survived==False):

[~/.local/lib/python3.7/site-packages/taurunner/Casino.py](https://vscode-remote+ssh-002dremote-002blogin-002erc-002efas-002eharvard-002eedu.vscode-resource.vscode-cdn.net/n/home11/tbertolez/IceCubeZoomEeV/~/.local/lib/python3.7/site-packages/taurunner/Casino.py) in Propagate(particle, track, body, condition)
     20     particle : taurunner Particle object after propagation
     21     '''
---> 22     total_column_depth = track.total_column_depth(body)
     23     total_distance     = track.x_to_d(1.-particle.position)*body.radius/units.km
     24     #keep iterating until final column depth is reached or a charged lepton is made

[~/.local/lib/python3.7/site-packages/taurunner/track/track.py](https://vscode-remote+ssh-002dremote-002blogin-002erc-002efas-002eharvard-002eedu.vscode-resource.vscode-cdn.net/n/home11/tbertolez/IceCubeZoomEeV/~/.local/lib/python3.7/site-packages/taurunner/track/track.py) in total_column_depth(self, body, safe_mode)
     51         hash_s = get_hash(self, body)
     52         if hash_s not in self._column_depth_lookup.keys():
---> 53             set_spline(self, body)
     54         return self._column_depth_lookup[hash_s][0]
     55 

[~/.local/lib/python3.7/site-packages/taurunner/track/utils/spline_column_depth.py](https://vscode-remote+ssh-002dremote-002blogin-002erc-002efas-002eharvard-002eedu.vscode-resource.vscode-cdn.net/n/home11/tbertolez/IceCubeZoomEeV/~/.local/lib/python3.7/site-packages/taurunner/track/utils/spline_column_depth.py) in set_spline(track, body, npad)
    150 
    151     # Construct the X2x on the fly
--> 152     X2x = construct_X2x(x2X)
    153 
    154     track._column_depth_lookup[hash_s] = (tot_X, x2X, X2x)

[~/.local/lib/python3.7/site-packages/taurunner/track/utils/spline_column_depth.py](https://vscode-remote+ssh-002dremote-002blogin-002erc-002efas-002eharvard-002eedu.vscode-resource.vscode-cdn.net/n/home11/tbertolez/IceCubeZoomEeV/~/.local/lib/python3.7/site-packages/taurunner/track/utils/spline_column_depth.py) in construct_X2x(x2X)
    105         cds        = [float(x2X(x)) for x in xx]
    106         cds[0]     = 0. # Sometimes spline support can be imperfect
--> 107         X2x_tck    = splrep(cds, xx)
    108         def X2x(X):
    109             return float(splev(X, X2x_tck))

/cvmfs/icecube.opensciencegrid.org/py3-v4.1.1/RHEL_8_x86_64/lib/python3.7/site-packages/scipy/interpolate/fitpack.py in splrep(x, y, w, xb, xe, k, task, s, t, full_output, per, quiet)
    287 
    288     """
--> 289     res = _impl.splrep(x, y, w, xb, xe, k, task, s, t, full_output, per, quiet)
    290     return res
    291 

/cvmfs/icecube.opensciencegrid.org/py3-v4.1.1/RHEL_8_x86_64/lib/python3.7/site-packages/scipy/interpolate/_fitpack_impl.py in splrep(x, y, w, xb, xe, k, task, s, t, full_output, per, quiet)
    513         else:
    514             try:
--> 515                 raise _iermess[ier][1](_iermess[ier][0])
    516             except KeyError:
    517                 raise _iermess['unknown'][1](_iermess['unknown'][0])

ValueError: Error on input data

If I write th_max = 179.9, then I recover the same error message as William states. If I set th_max = 90, everything runs smoothly! Hope this helps :)