#include <math.h>
#include "calc_mag.h"

double a0_dist[]={0,5,10,20,35,55,73,96,220,334,445,600,750,1010,1340,1600,1750,1960,2230};
double a0_val[]={1.4,1.4,1.5,1.7,2.3,2.7,2.85,3.0,3.65,4.2,4.6,5.0,5.3,5.7,6,6.2,6.3,6.4,6.5};

double log_a0(double dist) {
  int i;
  int found=-1;
  int prev;
  double ret_val;

  for (i=0;i<19 &&found ==-1 ;i++) {
    if (a0_dist[i]>=dist) {
      found=i;
    }
  }
  prev=found-1;

  if (a0_dist[found]==dist) {
    ret_val= a0_val[found];
  } else if (found==-1) {
    ret_val= a0_val[19];
  } else {
    ret_val=-(a0_val[prev]+(a0_val[found]-a0_val[prev])*(dist-a0_dist[prev])/(a0_dist[found]-a0_dist[prev]));
  }
  return ret_val;
}

double calc_ml (double amplitude, Station_Struct sta_info, double lat, double lon) {
  double delta=calc_dist(sta_info.lat,sta_info.lon,lat,lon)*111.195;
  double ret_val=0;
  
  ret_val=log10(amplitude)-log_a0(delta)+sta_info.ml;

  return ret_val;
}