#include "geodesic.h"
#include "observer_track.h"

#include <math.h>
#include <stdio.h>
#include <string.h>

static int nearly_equal(double a, double b) { return fabs(a - b) < 1e-12; }

/* Every enumerator must have a stable label, appended detail codes must be
 * distinct and valid, the coarse category helper must stay inside the coarse
 * range, and both out-of-range directions (negative and >= COUNT) must be
 * rejected as UNKNOWN/sentinel.  The frozen wire ids 0..9 are pinned. */
static int check_reason_names(void) {
  static const char *const frozen[] = {
      "NONE",                "REDSHIFT_LIMIT",   "BUDGET_EXHAUSTED",
      "TIME_RANGE_EXHAUSTED", "OUT_OF_DOMAIN",   "INVALID_METRIC",
      "INTEGRATION_ERROR",   "UNSUPPORTED",      "PROTOCOL_ERROR",
      "IO_ERROR"};
  int failed = 0;
  for (size_t i = 0; i < sizeof frozen / sizeof frozen[0]; ++i) {
    const RayReason r = (RayReason)i;
    if (!ray_reason_valid(r) || strcmp(ray_reason_name(r), frozen[i]) != 0 ||
        ray_reason_category(r) != r) {
      fprintf(stderr, "frozen reason id %zu is not stable\n", i);
      failed = 1;
    }
  }
  if ((unsigned)RAY_REASON_IO_ERROR + 1u >= (unsigned)RAY_REASON_COUNT) {
    fputs("no appended detail reasons\n", stderr);
    failed = 1;
  }
  for (unsigned i = 0; i < (unsigned)RAY_REASON_COUNT; ++i) {
    const RayReason r = (RayReason)i;
    const char *name = ray_reason_name(r);
    if (!ray_reason_valid(r) || name == NULL || name[0] == '\0' ||
        strcmp(name, "UNKNOWN") == 0) {
      fprintf(stderr, "reason %u has no stable label\n", i);
      failed = 1;
    }
    if (!ray_reason_valid(ray_reason_category(r))) {
      fprintf(stderr, "reason %u has no valid category\n", i);
      failed = 1;
    }
  }
  if (ray_reason_category(RAY_REASON_WORLDTUBE_SAMPLE_FAILED) !=
          RAY_REASON_PROTOCOL_ERROR ||
      ray_reason_category(RAY_REASON_OUTSIDE_WORLDTUBE) !=
          RAY_REASON_PROTOCOL_ERROR ||
      ray_reason_category(RAY_REASON_ESCAPE_LOCALIZATION_FAILED) !=
          RAY_REASON_PROTOCOL_ERROR ||
      ray_reason_category(RAY_REASON_REJECTION_LIMIT) !=
          RAY_REASON_INTEGRATION_ERROR ||
      ray_reason_category(RAY_REASON_INVALID_ESCAPE_DIRECTION) !=
          RAY_REASON_INTEGRATION_ERROR ||
      ray_reason_category(RAY_REASON_ENTRY_UNCONFIRMED) !=
          RAY_REASON_INTEGRATION_ERROR ||
      ray_reason_category(RAY_REASON_SLAB_LOAD_FAILED) != RAY_REASON_IO_ERROR) {
    fputs("detail reasons map to the wrong coarse category\n", stderr);
    failed = 1;
  }
  if (ray_reason_valid((RayReason)RAY_REASON_COUNT) ||
      ray_reason_valid((RayReason)-1)) {
    fputs("out-of-range reason accepted as valid\n", stderr);
    failed = 1;
  }
  if (strcmp(ray_reason_name((RayReason)RAY_REASON_COUNT), "UNKNOWN") != 0) {
    fputs("sentinel reason name is not UNKNOWN\n", stderr);
    failed = 1;
  }
  const RayReason unknown = (RayReason)((unsigned)RAY_REASON_COUNT + 7u);
  if (strcmp(ray_reason_name(unknown), "UNKNOWN") != 0 ||
      ray_reason_valid(unknown) ||
      ray_reason_category(unknown) != (RayReason)RAY_REASON_COUNT) {
    fputs("unknown reason is not the UNKNOWN sentinel\n", stderr);
    failed = 1;
  }
  return failed;
}

/* A NULL, non-finite or non-unit direction must be rejected as INVALID_ARGUMENT
 * without dereferencing the direction or running any trace. */
static int check_invalid_directions(const SpacetimeSource *source,
                                    const ObserverState *observer) {
  const GeodesicTraceConfig config = {.coordinate_time_step = 0.25,
                                      .max_steps = 100};
  int failed = 0;
  const RayEndpoint null_dir =
      geodesic_trace_past(source, observer, NULL, &config);
  if (null_dir.outcome != RAY_OUTCOME_INCOMPLETE ||
      null_dir.reason != RAY_REASON_INVALID_ARGUMENT) {
    fputs("NULL direction is not INVALID_ARGUMENT\n", stderr);
    failed = 1;
  }
  const double nan_dir[3] = {NAN, 0.0, 0.0};
  const RayEndpoint nan =
      geodesic_trace_past(source, observer, nan_dir, &config);
  if (nan.outcome != RAY_OUTCOME_INCOMPLETE ||
      nan.reason != RAY_REASON_INVALID_ARGUMENT) {
    fputs("NaN direction is not INVALID_ARGUMENT\n", stderr);
    failed = 1;
  }
  const double nonunit_dir[3] = {2.0, 0.0, 0.0};
  const RayEndpoint nonunit =
      geodesic_trace_past(source, observer, nonunit_dir, &config);
  if (nonunit.outcome != RAY_OUTCOME_INCOMPLETE ||
      nonunit.reason != RAY_REASON_INVALID_ARGUMENT) {
    fputs("non-unit direction is not INVALID_ARGUMENT\n", stderr);
    failed = 1;
  }
  return failed;
}

static int check_ray(const SpacetimeSource *source,
                     const ObserverState *observer,
                     const double local_direction[3],
                     const double expected[3]) {  const GeodesicTraceConfig config = {.coordinate_time_step = 0.25,
                                      .max_steps = 100};
  RayEndpoint ray =
      geodesic_trace_past(source, observer, local_direction, &config);
  if (ray.outcome != RAY_OUTCOME_ESCAPED ||
      !nearly_equal(ray.frequency_ratio, 1.0) ||
      !nearly_equal(ray.n_infinity[0], expected[0]) ||
      !nearly_equal(ray.n_infinity[1], expected[1]) ||
      !nearly_equal(ray.n_infinity[2], expected[2])) {
    fprintf(stderr, "flat-space geodesic regression failed\n");
    return 1;
  }
  return 0;
}

int main(void) {
  SpacetimeSource source = {0};
  MetricSlab *slab = NULL;
  MetricData metric;
  const ObserverState observer = observer_fixed_at_origin();
  if (spacetime_create_minkowski(&source, 10.0) ||
      spacetime_load_slab(&source, 0.0, -1.0, &slab) ||
      spacetime_slab_eval(slab, -0.5, (double[]){0.0, 0.0, 0.0}, &metric) ||
      metric.alpha != 1.0 ||
      !spacetime_slab_eval(slab, 0.25, (double[]){0.0, 0.0, 0.0}, &metric))
    return 1;
  int result = check_reason_names() ||
               check_invalid_directions(&source, &observer) ||
               check_ray(&source, &observer, (double[]){1.0, 0.0, 0.0},
                         (double[]){0.0, 0.0, -1.0}) ||
               check_ray(&source, &observer, (double[]){0.0, 0.0, 1.0},
                         (double[]){1.0, 0.0, 0.0});
  const ObserverState look_at_ra_zero =
      observer_fixed_at_origin_look_at(0.0, 0.0);
  result = result || check_ray(&source, &look_at_ra_zero,
                               (double[]){1.0, 0.0, 0.0},
                               (double[]){1.0, 0.0, 0.0});
  /* Standard ICRS has +Z at the north celestial pole and +Y at increasing
   * RA.  A right-handed north-up camera consequently has west to its right. */
  result = result || check_ray(&source, &look_at_ra_zero,
                               (double[]){0.0, 1.0, 0.0},
                               (double[]){0.0, 0.0, 1.0}) ||
           check_ray(&source, &look_at_ra_zero,
                     (double[]){0.0, 0.0, 1.0},
                     (double[]){0.0, -1.0, 0.0});
  ObserverTrack accelerated = {0};
  ObserverState final_observer;
  if (observer_track_generate_minkowski_acceleration(&accelerated, 1.52, 2.0,
                                                      1.0 / 30.0) ||
      observer_track_interpolate(&accelerated, 2.0, &final_observer, NULL)) {
    result = 1;
  } else {
    const RayEndpoint forward = geodesic_trace_past(
        &source, &final_observer, (double[]){1.0, 0.0, 0.0},
        &(GeodesicTraceConfig){.coordinate_time_step = 0.25, .max_steps = 100});
    const double expected_g = sqrt(1.0 + 3.04 * 3.04) + 3.04;
    if (forward.outcome != RAY_OUTCOME_ESCAPED ||
        !nearly_equal(forward.frequency_ratio, expected_g)) {
      fputs("accelerated-observer Doppler regression failed\n", stderr);
      result = 1;
    }
  }
  observer_track_destroy(&accelerated);
  spacetime_free_slab(slab);
  spacetime_destroy(&source);
  return result;
}