def METHOD_NAME(zenith, azimuth, positions, n=1.000293):
    """
    Calculate the time delay between given positions for an arrival direction
    Parameters
    ----------
    zenith: float [rad]
        Zenith angle in convention up = 0
    azimuth: float [rad]
        Azimuth angle in convention East = 0, counter-clock-wise
    positions: array[N x 3]
        Positions on ground
    n: float (default: 1.000293)
        Index of reflection of propagation medium. By default, air is assumed
    """
    shower_axis = np.array([np.sin(zenith) * np.cos(azimuth), np.sin(zenith) * np.sin(azimuth), np.cos(zenith)])
    if positions.ndim == 1:
        times = -(1 / (constants.c / n)) * np.dot(shower_axis, positions) * units.s
    else:
        times = np.zeros(positions.shape[0])
        for i in range(positions.shape[0]):
            times[i] = -(1 / (constants.c / n)) * np.dot(shower_axis, positions[i, :]) * units.s
    return times