/* Copyright (c) 2019, Anthony Latorre <tlatorre at uchicago>
 *
 * This program is free software: you can redistribute it and/or modify it
 * under the terms of the GNU General Public License as published by the Free
 * Software Foundation, either version 3 of the License, or (at your option)
 * any later version.

 * This program is distributed in the hope that it will be useful, but WITHOUT
 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for
 * more details.

 * You should have received a copy of the GNU General Public License along with
 * this program. If not, see <https://www.gnu.org/licenses/>.
 */

#ifndef PATH_H
#define PATH_H

#include <gsl/gsl_spline.h>
#include <stddef.h> /* for size_t */

/* Minimum value for the scattering RMS `theta0`.
 *
 * This is a bit of a hack to correct for the fact that we assume the
 * probability of a photon hitting a PMT is uniform across the face of a PMT.
 * By introducing a minimum value for the scattering RMS we broaden the angular
 * distribution such that it effectively averages across the face of a PMT.
 *
 * Initially had this at 0.01 since that roughly corresponded to the angular
 * width of the PMT from the center to the PSUP. Now, since we calculate the
 * angular width for every point in get_theta0_min() this isn't necessary. We
 * still keep it though so that at the *very* beginning of the track theta0 is
 * not zero. */
#define MIN_THETA0 1e-5

typedef double getKineticEnergyFunc(double x, void *params);

typedef struct path {
    double pos[3];
    double dir[3];
    double theta0;

    size_t n;
    size_t len;

    double *s;
    double *x;
    double *y;
    double *z;
    double *beta;
    double *t;
    double *dx;
    double *dy;
    double *dz;
} path;

double path_get_coefficient(unsigned int k, double *s, double *x, double theta0, size_t n);
path *path_init(double *pos, double *dir, double T0, double range, size_t len, double theta0, getKineticEnergyFunc *fun, void *params, double *z1, double *z2, size_t n, double m);
void path_eval(path *p, double s, double *pos, double *dir, double *beta, double *t, double *theta0);
void path_free(path *p);

#endif