#include <math.h>
#include <seiscomp3/seismology/magnitudes.h>

static float __qmb[17][108] = {
    { 
0.0, 0.0, 0.0, 5.75, 5.88, 5.96, 6.01, 6.10, 6.14, 6.17, 6.17, 6.16, 6.15, 6.15, 6.10, 6.06, 6.05, 6.04, 6.05, 6.08, 6.11, 6.16, 6.23, 6.30, 6.39, 6.48, 6.55, 6.61, 6.64, 6.66, 6.67, 6.67, 6.67, 6.67, 6.66, 6.66, 6.65, 6.66, 6.66, 6.67, 6.67, 6.69, 6.69, 6.70, 6.71, 6.72, 6.72, 6.73, 6.74, 6.75, 6.76, 6.77, 6.77, 6.78, 6.79, 6.79, 6.80, 6.81, 6.82, 6.82, 6.83, 6.83, 6.84, 6.84, 6.84, 6.85, 6.85, 6.85, 6.86, 6.87, 6.87, 6.87, 6.87, 6.88, 6.88, 6.88, 6.89, 6.89, 6.90, 6.89, 6.89, 6.89, 6.90, 6.92, 6.94, 6.97, 7.01, 7.05, 7.09, 7.15, 7.20, 7.25, 7.31, 7.37, 7.42, 7.46, 7.50, 7.55, 7.59, 7.64, 7.68, 7.73, 7.75, 7.78, 7.78, 7.78, 7.78, 7.78
    },
    { 
0.0, 0.0, 0.0, 5.77, 5.85, 5.92, 5.99, 6.09, 6.16, 6.21, 6.23, 6.22, 6.20, 6.19, 6.15, 6.12, 6.10, 6.11, 6.12, 6.14, 6.17, 6.21, 6.26, 6.33, 6.41, 6.50, 6.57, 6.63, 6.65, 6.68, 6.69, 6.69, 6.69, 6.69, 6.69, 6.69, 6.69, 6.69, 6.68, 6.68, 6.68, 6.68, 6.69, 6.70, 6.71, 6.72, 6.74, 6.74, 6.74, 6.75, 6.75, 6.76, 6.77, 6.77, 6.78, 6.78, 6.78, 6.78, 6.79, 6.79, 6.80, 6.81, 6.81, 6.81, 6.82, 6.82, 6.83, 6.83, 6.83, 6.83, 6.83, 6.83, 6.83, 6.84, 6.84, 6.84, 6.85, 6.85, 6.85, 6.86, 6.86, 6.86, 6.86, 6.88, 6.90, 6.94, 6.98, 7.02, 7.06, 7.12, 7.16, 7.21, 7.26, 7.31, 7.36, 7.41, 7.45, 7.50, 7.54, 7.58, 7.62, 7.65, 7.67, 7.68, 7.68, 7.68, 7.68, 7.68
    },
    { 
0.0, 0.0, 0.0, 5.70, 5.86, 5.92, 5.96, 6.03, 6.10, 6.12, 6.15, 6.17, 6.16, 6.14, 6.13, 6.10, 6.08, 6.10, 6.12, 6.14, 6.16, 6.21, 6.24, 6.32, 6.40, 6.49, 6.54, 6.60, 6.63, 6.66, 6.69, 6.70, 6.69, 6.68, 6.68, 6.65, 6.65, 6.66, 6.66, 6.65, 6.66, 6.65, 6.65, 6.66, 6.68, 6.69, 6.71, 6.73, 6.74, 6.75, 6.76, 6.75, 6.74, 6.74, 6.75, 6.75, 6.76, 6.77, 6.77, 6.76, 6.76, 6.76, 6.76, 6.76, 6.78, 6.77, 6.79, 6.78, 6.78, 6.78, 6.79, 6.78, 6.78, 6.79, 6.79, 6.78, 6.78, 6.78, 6.78, 6.79, 6.79, 6.79, 6.80, 6.81, 6.83, 6.87, 6.92, 6.97, 7.02, 7.07, 7.12, 7.16, 7.21, 7.27, 7.33, 7.38, 7.41, 7.45, 7.49, 7.51, 7.54, 7.56, 7.58, 7.60, 7.60, 7.60, 7.60, 7.60
    },
    { 
0.0, 0.0, 0.0, 5.73, 5.80, 5.82, 5.88, 6.01, 6.10, 6.14, 6.15, 6.14, 6.14, 6.09, 6.04, 6.03, 6.03, 6.03, 6.07, 6.11, 6.13, 6.18, 6.22, 6.27, 6.37, 6.46, 6.53, 6.59, 6.62, 6.64, 6.66, 6.67, 6.67, 6.66, 6.65, 6.63, 6.64, 6.63, 6.63, 6.64, 6.63, 6.63, 6.64, 6.65, 6.66, 6.66, 6.67, 6.69, 6.72, 6.74, 6.75, 6.74, 6.74, 6.72, 6.72, 6.73, 6.72, 6.72, 6.73, 6.72, 6.73, 6.73, 6.74, 6.74, 6.74, 6.73, 6.74, 6.74, 6.74, 6.74, 6.75, 6.75, 6.76, 6.77, 6.78, 6.78, 6.77, 6.76, 6.76, 6.76, 6.74, 6.75, 6.77, 6.79, 6.82, 6.87, 6.93, 6.99, 7.03, 7.08, 7.13, 7.17, 7.21, 7.28, 7.34, 7.39, 7.41, 7.47, 7.53, 7.57, 7.62, 7.63, 7.64, 7.64, 7.64, 7.64, 7.64, 7.64
    },
    { 
0.0, 0.0, 0.0, 5.67, 5.73, 5.73, 5.83, 5.96, 6.02, 6.10, 6.20, 6.18, 6.15, 6.10, 6.05, 6.03, 6.01, 6.03, 6.05, 6.09, 6.10, 6.14, 6.19, 6.25, 6.33, 6.43, 6.50, 6.56, 6.59, 6.61, 6.61, 6.61, 6.60, 6.59, 6.57, 6.57, 6.57, 6.57, 6.56, 6.59, 6.58, 6.59, 6.59, 6.60, 6.59, 6.60, 6.60, 6.60, 6.61, 6.63, 6.65, 6.66, 6.67, 6.67, 6.67, 6.67, 6.68, 6.68, 6.68, 6.66, 6.66, 6.66, 6.66, 6.66, 6.66, 6.66, 6.65, 6.65, 6.66, 6.66, 6.68, 6.68, 6.69, 6.69, 6.69, 6.70, 6.70, 6.70, 6.71, 6.71, 6.71, 6.72, 6.73, 6.75, 6.77, 6.82, 6.86, 6.92, 6.98, 7.03, 7.07, 7.12, 7.17, 7.22, 7.27, 7.32, 7.34, 7.38, 7.42, 7.47, 7.51, 7.54, 7.54, 7.56, 7.56, 7.56, 7.56, 7.56
    }, 
    { 
0.0, 0.0, 0.0, 5.64, 5.68, 5.70, 5.78, 5.90, 5.98, 6.03, 6.08, 6.03, 5.94, 5.88, 5.91, 5.92, 5.95, 6.00, 6.03, 6.04, 6.06, 6.09, 6.12, 6.18, 6.24, 6.33, 6.41, 6.47, 6.52, 6.51, 6.51, 6.49, 6.48, 6.45, 6.45, 6.45, 6.46, 6.46, 6.46, 6.47, 6.46, 6.46, 6.46, 6.47, 6.47, 6.50, 6.50, 6.52, 6.51, 6.52, 6.53, 6.53, 6.52, 6.54, 6.55, 6.55, 6.56, 6.56, 6.56, 6.55, 6.56, 6.56, 6.57, 6.56, 6.57, 6.57, 6.56, 6.56, 6.57, 6.57, 6.58, 6.59, 6.59, 6.59, 6.60, 6.60, 6.61, 6.62, 6.63, 6.63, 6.64, 6.65, 6.66, 6.68, 6.71, 6.75, 6.78, 6.84, 6.89, 6.94, 6.99, 7.04, 7.10, 7.14, 7.20, 7.27, 7.30, 7.33, 7.38, 7.42, 7.46, 7.49, 7.51, 7.53, 7.53, 7.53, 7.53, 7.53
    },
    { 
0.0, 0.0, 0.0, 5.49, 5.50, 5.56, 5.68, 5.76, 5.86, 5.85, 5.84, 5.73, 5.70, 5.68, 5.75, 5.77, 5.84, 5.89, 5.96, 5.94, 5.97, 6.01, 6.07, 6.13, 6.22, 6.29, 6.34, 6.38, 6.42, 6.41, 6.41, 6.40, 6.38, 6.35, 6.35, 6.36, 6.36, 6.37, 6.38, 6.39, 6.38, 6.37, 6.35, 6.34, 6.35, 6.37, 6.38, 6.41, 6.40, 6.41, 6.42, 6.42, 6.42, 6.44, 6.44, 6.42, 6.44, 6.44, 6.44, 6.46, 6.48, 6.47, 6.47, 6.45, 6.43, 6.43, 6.42, 6.43, 6.46, 6.47, 6.47, 6.49, 6.50, 6.51, 6.52, 6.53, 6.53, 6.54, 6.55, 6.54, 6.55, 6.55, 6.56, 6.58, 6.61, 6.64, 6.69, 6.73, 6.76, 6.81, 6.85, 6.89, 6.95, 7.03, 7.08, 7.17, 7.21, 7.22, 7.32, 7.35, 7.42, 7.45, 7.47, 7.47, 7.47, 7.47, 7.47, 7.47
    },
    { 
0.0, 0.0, 0.0, 5.49, 5.49, 5.57, 5.73, 5.82, 5.90, 5.85, 5.76, 5.66, 5.66, 5.64, 5.73, 5.76, 5.87, 5.88, 5.95, 5.94, 5.98, 6.02, 6.09, 6.15, 6.22, 6.26, 6.30, 6.35, 6.38, 6.38, 6.37, 6.35, 6.31, 6.28, 6.28, 6.29, 6.28, 6.30, 6.30, 6.28, 6.28, 6.28, 6.28, 6.30, 6.32, 6.32, 6.33, 6.34, 6.31, 6.33, 6.33, 6.35, 6.34, 6.38, 6.37, 6.38, 6.39, 6.38, 6.38, 6.40, 6.43, 6.43, 6.46, 6.43, 6.41, 6.39, 6.39, 6.39, 6.41, 6.44, 6.44, 6.45, 6.46, 6.47, 6.47, 6.50, 6.51, 6.53, 6.54, 6.56, 6.56, 6.56, 6.57, 6.59, 6.60, 6.63, 6.68, 6.73, 6.75, 6.80, 6.84, 6.89, 6.95, 7.03, 7.09, 7.17, 7.24, 7.27, 7.31, 7.35, 7.42, 7.43, 7.44, 7.47, 7.47, 7.47, 7.47, 7.47
    },
    { 
0.0, 0.0, 0.0, 5.52, 5.47, 5.55, 5.68, 5.76, 5.82, 5.84, 5.77, 5.70, 5.72, 5.72, 5.77, 5.77, 5.83, 5.84, 5.87, 5.88, 5.93, 5.99, 6.09, 6.20, 6.24, 6.27, 6.30, 6.31, 6.32, 6.32, 6.30, 6.28, 6.24, 6.23, 6.21, 6.21, 6.21, 6.23, 6.22, 6.22, 6.23, 6.22, 6.23, 6.24, 6.24, 6.25, 6.25, 6.23, 6.21, 6.24, 6.25, 6.28, 6.30, 6.34, 6.31, 6.34, 6.33, 6.31, 6.32, 6.35, 6.36, 6.38, 6.42, 6.41, 6.42, 6.42, 6.40, 6.39, 6.41, 6.43, 6.41, 6.42, 6.44, 6.44, 6.45, 6.48, 6.50, 6.52, 6.53, 6.54, 6.53, 6.53, 6.52, 6.55, 6.57, 6.61, 6.66, 6.70, 6.72, 6.76, 6.80, 6.84, 6.90, 6.97, 7.03, 7.11, 7.18, 7.22, 7.27, 7.30, 7.36, 7.38, 7.41, 7.42, 7.42, 7.42, 7.42, 7.42
    },
    { 
0.0, 0.0, 0.0, 5.58, 5.51, 5.55, 5.64, 5.74, 5.75, 5.77, 5.74, 5.72, 5.75, 5.78, 5.79, 5.79, 5.82, 5.82, 5.83, 5.87, 5.93, 5.99, 6.06, 6.15, 6.16, 6.18, 6.23, 6.27, 6.25, 6.27, 6.25, 6.21, 6.19, 6.21, 6.20, 6.20, 6.19, 6.18, 6.16, 6.17, 6.15, 6.16, 6.19, 6.20, 6.18, 6.21, 6.20, 6.20, 6.17, 6.20, 6.20, 6.24, 6.25, 6.31, 6.29, 6.32, 6.31, 6.29, 6.29, 6.30, 6.30, 6.35, 6.38, 6.38, 6.42, 6.42, 6.39, 6.40, 6.40, 6.40, 6.38, 6.37, 6.38, 6.38, 6.40, 6.43, 6.46, 6.49, 6.50, 6.52, 6.51, 6.51, 6.52, 6.56, 6.58, 6.61, 6.66, 6.71, 6.73, 6.77, 6.81, 6.85, 6.91, 6.97, 7.03, 7.09, 7.17, 7.22, 7.27, 7.30, 7.35, 7.36, 7.38, 7.38, 7.38, 7.38, 7.38, 7.38
    },
    { 
0.0, 0.0, 0.0, 5.44, 5.44, 5.46, 5.51, 5.61, 5.61, 5.71, 5.72, 5.74, 5.76, 5.79, 5.78, 5.78, 5.78, 5.80, 5.81, 5.86, 5.88, 5.99, 6.05, 6.15, 6.15, 6.16, 6.18, 6.21, 6.18, 6.19, 6.19, 6.17, 6.14, 6.17, 6.15, 6.14, 6.14, 6.14, 6.12, 6.14, 6.13, 6.13, 6.16, 6.18, 6.15, 6.16, 6.16, 6.17, 6.14, 6.18, 6.18, 6.21, 6.22, 6.28, 6.27, 6.28, 6.28, 6.26, 6.27, 6.26, 6.27, 6.30, 6.32, 6.32, 6.38, 6.38, 6.36, 6.38, 6.38, 6.37, 6.36, 6.37, 6.37, 6.39, 6.40, 6.42, 6.44, 6.46, 6.46, 6.47, 6.48, 6.49, 6.49, 6.53, 6.57, 6.62, 6.66, 6.71, 6.74, 6.77, 6.81, 6.86, 6.92, 6.99, 7.06, 7.12, 7.20, 7.24, 7.29, 7.32, 7.36, 7.37, 7.41, 7.39, 7.39, 7.39, 7.39, 7.39
    },
    { 
0.0, 0.0, 0.0, 5.38, 5.43, 5.49, 5.52, 5.57, 5.59, 5.66, 5.73, 5.76, 5.77, 5.78, 5.79, 5.81, 5.82, 5.84, 5.85, 5.89, 5.95, 6.02, 6.06, 6.11, 6.15, 6.18, 6.18, 6.17, 6.17, 6.18, 6.16, 6.16, 6.15, 6.14, 6.15, 6.14, 6.14, 6.14, 6.14, 6.14, 6.14, 6.16, 6.15, 6.14, 6.14, 6.14, 6.14, 6.16, 6.17, 6.18, 6.19, 6.21, 6.23, 6.23, 6.24, 6.25, 6.26, 6.26, 6.26, 6.27, 6.27, 6.29, 6.31, 6.31, 6.32, 6.33, 6.34, 6.35, 6.34, 6.35, 6.35, 6.36, 6.37, 6.38, 6.39, 6.41, 6.41, 6.42, 6.44, 6.44, 6.46, 6.47, 6.49, 6.52, 6.54, 6.58, 6.62, 6.66, 6.69, 6.73, 6.78, 6.81, 6.87, 6.92, 6.99, 7.04, 7.10, 7.14, 7.19, 7.22, 7.25, 7.27, 7.28, 7.32, 7.32, 7.32, 7.32, 7.32
    },
    { 
0.0, 0.0, 0.0, 5.37, 5.40, 5.44, 5.48, 5.54, 5.57, 5.65, 5.71, 5.75, 5.75, 5.78, 5.78, 5.80, 5.82, 5.85, 5.85, 5.90, 5.94, 6.00, 6.04, 6.08, 6.12, 6.14, 6.14, 6.14, 6.14, 6.15, 6.15, 6.14, 6.13, 6.13, 6.14, 6.14, 6.12, 6.11, 6.11, 6.11, 6.11, 6.12, 6.11, 6.11, 6.11, 6.12, 6.13, 6.15, 6.16, 6.17, 6.18, 6.20, 6.21, 6.21, 6.22, 6.23, 6.23, 6.23, 6.23, 6.25, 6.25, 6.26, 6.30, 6.29, 6.31, 6.32, 6.34, 6.34, 6.34, 6.34, 6.34, 6.35, 6.36, 6.37, 6.37, 6.38, 6.38, 6.39, 6.41, 6.42, 6.44, 6.45, 6.47, 6.50, 6.52, 6.57, 6.61, 6.65, 6.68, 6.72, 6.77, 6.81, 6.86, 6.93, 6.99, 7.04, 7.10, 7.13, 7.18, 7.20, 7.23, 7.25, 7.28, 7.34, 7.34, 7.34, 7.34, 7.34
    },
    { 
0.0, 0.0, 0.0, 5.42, 5.44, 5.48, 5.50, 5.55, 5.57, 5.63, 5.66, 5.70, 5.72, 5.75, 5.78, 5.80, 5.84, 5.88, 5.91, 5.95, 5.98, 6.05, 6.10, 6.13, 6.15, 6.16, 6.15, 6.14, 6.12, 6.13, 6.12, 6.12, 6.11, 6.11, 6.11, 6.10, 6.10, 6.10, 6.10, 6.10, 6.09, 6.10, 6.09, 6.09, 6.10, 6.11, 6.13, 6.15, 6.17, 6.19, 6.20, 6.22, 6.22, 6.21, 6.21, 6.21, 6.22, 6.21, 6.21, 6.23, 6.23, 6.25, 6.28, 6.28, 6.29, 6.30, 6.31, 6.30, 6.30, 6.31, 6.30, 6.31, 6.32, 6.33, 6.34, 6.35, 6.37, 6.38, 6.40, 6.41, 6.43, 6.45, 6.47, 6.50, 6.53, 6.57, 6.61, 6.65, 6.68, 6.73, 6.77, 6.81, 6.86, 6.92, 6.96, 7.01, 7.07, 7.11, 7.15, 7.18, 7.20, 7.23, 7.25, 7.29, 7.29, 7.29, 7.29, 7.29
    },
    { 
0.0, 0.0, 0.0, 5.43, 5.44, 5.48, 5.50, 5.53, 5.53, 5.59, 5.61, 5.65, 5.69, 5.74, 5.78, 5.82, 5.88, 5.94, 6.01, 6.02, 6.08, 6.15, 6.20, 6.22, 6.21, 6.21, 6.20, 6.16, 6.13, 6.13, 6.11, 6.11, 6.12, 6.12, 6.11, 6.10, 6.10, 6.11, 6.11, 6.10, 6.11, 6.09, 6.10, 6.10, 6.10, 6.12, 6.14, 6.16, 6.17, 6.18, 6.19, 6.21, 6.21, 6.20, 6.19, 6.19, 6.19, 6.19, 6.20, 6.21, 6.21, 6.22, 6.25, 6.26, 6.27, 6.27, 6.28, 6.26, 6.27, 6.27, 6.26, 6.26, 6.28, 6.29, 6.30, 6.31, 6.34, 6.36, 6.37, 6.39, 6.42, 6.45, 6.48, 6.51, 6.53, 6.56, 6.60, 6.64, 6.68, 6.72, 6.76, 6.80, 6.85, 6.91, 6.96, 7.01, 7.07, 7.11, 7.14, 7.18, 7.19, 7.22, 7.25, 7.28, 7.28, 7.28, 7.28, 7.28
    },
    { 
0.0, 0.0, 0.0, 5.41, 5.31, 5.38, 5.37, 5.37, 5.43, 5.45, 5.51, 5.58, 5.64, 5.69, 5.76, 5.85, 5.89, 5.98, 6.10, 6.09, 6.19, 6.24, 6.23, 6.21, 6.20, 6.20, 6.16, 6.11, 6.10, 6.07, 6.03, 6.05, 6.10, 6.10, 6.10, 6.13, 6.14, 6.15, 6.13, 6.14, 6.15, 6.14, 6.15, 6.14, 6.13, 6.15, 6.16, 6.18, 6.17, 6.18, 6.18, 6.19, 6.20, 6.18, 6.19, 6.19, 6.19, 6.20, 6.20, 6.21, 6.21, 6.21, 6.25, 6.24, 6.26, 6.27, 6.28, 6.28, 6.28, 6.27, 6.26, 6.25, 6.27, 6.28, 6.29, 6.29, 6.30, 6.32, 6.33, 6.35, 6.39, 6.43, 6.47, 6.49, 6.52, 6.55, 6.58, 6.64, 6.69, 6.73, 6.77, 6.82, 6.88, 6.92, 6.96, 7.03, 7.08, 7.10, 7.12, 7.16, 7.20, 7.21, 7.25, 7.28, 7.28, 7.28, 7.28, 7.28
    },
    { 
0.0, 0.0, 0.0, 4.89, 5.10, 5.29, 5.29, 5.27, 5.27, 5.44, 5.64, 5.67, 5.67, 5.79, 5.79, 5.89, 6.17, 6.05, 6.14, 6.16, 6.13, 6.16, 6.11, 6.08, 6.09, 6.07, 6.04, 6.00, 6.03, 6.01, 6.03, 6.05, 6.05, 6.08, 6.12, 6.15, 6.17, 6.21, 6.23, 6.21, 6.22, 6.19, 6.17, 6.17, 6.16, 6.16, 6.12, 6.10, 6.09, 6.10, 6.12, 6.11, 6.10, 6.11, 6.11, 6.11, 6.09, 6.10, 6.19, 6.20, 6.21, 6.20, 6.23, 6.20, 6.21, 6.19, 6.23, 6.22, 6.22, 6.20, 6.20, 6.23, 6.24, 6.27, 6.28, 6.33, 6.35, 6.39, 6.42, 6.47, 6.53, 6.54, 6.58, 6.58, 6.60, 6.62, 6.64, 6.69, 6.70, 6.75, 6.74, 6.76, 6.83, 6.88, 6.96, 6.96, 7.00, 7.03, 7.07, 7.13, 7.18, 7.18, 7.21, 7.24, 7.24, 7.24, 7.24, 7.24
    }
};

#define max(x,y) ((x)>(y)?(x):(y))
#define min(x,y) ((x)<(y)?(x):(y))

static double
bmagnz(double amplitude, double period, double delta, double depth, int *err)
{
	int j, k;
	double mb, s1, s2, q1, Q;

	*err = 0;

	if (depth < 0. || depth > 700. || delta < 5 || delta > 105) {
		*err = -1;
		return 0.;
	}

	if (depth < 100) {
		k  = (int) (depth/25.) + 2;
		s1 = 0.04*(depth-25.*(k-2));
	}
	else {
		k = min((int)(depth/50.)+4, 17);
		s1 = .02*(depth-50.*(k-4));
	}
	j  = max(min((int)(delta), 108), 2);
	s2 = delta-j;
	q1 = __qmb[k-2][j-2] + s1*(__qmb[k-1][j-2]-__qmb[k-2][j-2]);
	Q  = q1+s2*((__qmb[k-2][j-1]+s1*(__qmb[k-1][j-1]-__qmb[k-2][j-1]))-q1);

	mb = log10(amplitude/period) + Q;

	return mb;
}

namespace Seiscomp {
namespace Magnitudes {

bool
compute_mb(
	double amplitude, 
	double period,    
	double delta,     
	double depth,     
	double *mag)
{
	int err;
	double m = bmagnz(amplitude, period, delta, depth, &err);
	if (err) return false;

	*mag = m;

	return true;
}

bool
compute_mb_fromVelocity(
	double amplitude, 
	double delta,     
	double depth,     
	double *mag)
{
	int err;
	double m = bmagnz(amplitude, 2*M_PI, delta, depth, &err);
	if (err) return false;

	*mag = m;

	return true;
}

} 
}
/***************************************************************************
 *   Copyright (C) by GFZ Potsdam                                          *
 *                                                                         *
 *   You can redistribute and/or modify this program under the             *
 *   terms of the SeisComP Public License.                                 *
 *                                                                         *
 *   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         *
 *   SeisComP Public License for more details.                             *
 ***************************************************************************/



#include <math.h>
#include <seiscomp3/seismology/magnitudes.h>

/* See http://jclahr.com/science/software/magnitude/mb/index.html
 *
 * These are the original Gutenberg and Richter (1956) calibration functions
 * as they are still (2007) used at the NEIC
 *
 * The Q values for a distance of 5 degrees are in the 4th column, i.e.
 * at __qmb[x][3]
 */

/*
static float __qmb[17][108] = {
    { 5.6, 5.8, 6.1, 6.4, 6.5, 7.0, 7.0, 7.2, 7.3, 7.2, 7.1, 7.0, 6.6, 6.3, 5.9, 5.9, 5.9, 6.0, 6.1, 6.1, 6.2, 6.3, 6.4, 6.5, 6.5, 6.5, 6.6, 6.6, 6.6, 6.7, 6.7, 6.7, 6.7, 6.6, 6.6, 6.5, 6.5, 6.4, 6.4, 6.5, 6.5, 6.5, 6.6, 6.7, 6.8, 6.9, 6.9, 6.8, 6.7, 6.7, 6.7, 6.7, 6.8, 6.8, 6.8, 6.8, 6.8, 6.9, 6.9, 6.9, 7.0, 7.0, 7.0, 7.0, 7.0, 7.0, 7.0, 7.0, 6.9, 6.9, 6.9, 6.9, 6.8, 6.8, 6.9, 6.9, 6.9, 6.8, 6.7, 6.8, 6.9, 7.0, 7.0, 7.0, 6.9, 7.0, 7.1, 7.0, 7.0, 7.1, 7.1, 7.2, 7.1, 7.2, 7.3, 7.4, 7.5, 7.5, 7.3, 7.4, 7.4, 7.5, 7.6, 7.7, 7.8, 7.8, 7.9, 8.0 },
    { 0.0, 0.0, 0.0, 6.3, 6.5, 6.8, 7.0, 7.0, 7.1, 7.0, 7.0, 6.9, 6.5, 6.1, 5.9, 5.9, 5.9, 6.0, 6.1, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.4, 6.5, 6.5, 6.6, 6.6, 6.7, 6.7, 6.7, 6.7, 6.7, 6.6, 6.6, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.7, 6.7, 6.8, 6.8, 6.8, 6.8, 6.7, 6.7, 6.7, 6.8, 6.8, 6.8, 6.8, 6.8, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.9, 7.0, 7.0, 7.0, 7.0, 7.1, 7.1, 7.0, 7.1, 7.2, 7.2, 7.2, 7.2, 7.2, 7.3, 7.3, 7.3, 7.3, 7.3, 7.4, 7.5, 7.6, 7.7, 7.8, 7.8, 7.9, 8.0 },
    { 0.0, 0.0, 0.0, 6.0, 6.2, 6.5, 6.6, 6.8, 6.9, 6.9, 6.8, 6.7, 6.5, 6.0, 5.9, 5.9, 5.9, 6.0, 6.1, 6.1, 6.2, 6.2, 6.3, 6.3, 6.3, 6.4, 6.4, 6.4, 6.5, 6.5, 6.6, 6.6, 6.7, 6.7, 6.7, 6.7, 6.7, 6.6, 6.6, 6.5, 6.5, 6.5, 6.5, 6.6, 6.7, 6.7, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.8, 6.8, 6.8, 6.8, 6.8, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.8, 6.8, 6.7, 6.7, 6.7, 6.8, 6.8, 6.8, 6.9, 7.0, 7.0, 7.0, 7.1, 7.1, 7.2, 7.2, 7.2, 7.2, 7.2, 7.3, 7.3, 7.3, 7.3, 7.3, 7.4, 7.4, 7.5, 7.6, 7.7, 7.8, 7.9, 7.9, 8.0 },
    { 0.0, 0.0, 0.0, 6.0, 6.1, 6.4, 6.5, 6.6, 6.8, 6.7, 6.7, 6.5, 6.1, 6.0, 5.9, 6.0, 6.0, 6.0, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.5, 6.5, 6.6, 6.6, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.6, 6.6, 6.6, 6.6, 6.6, 6.7, 6.7, 6.7, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.9, 6.9, 6.9, 6.9, 6.9, 6.9, 6.8, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.8, 6.8, 6.8, 6.9, 6.9, 7.0, 7.0, 7.1, 7.2, 7.2, 7.2, 7.2, 7.3, 7.3, 7.3, 7.3, 7.4, 7.4, 7.5, 7.5, 7.7, 7.8, 7.8, 7.9, 7.9, 8.0 },
    { 0.0, 0.0, 0.0, 6.0, 6.0, 6.3, 6.4, 6.6, 6.7, 6.6, 6.5, 6.4, 6.0, 6.0, 6.0, 6.0, 6.0, 6.1, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.5, 6.5, 6.5, 6.6, 6.6, 6.7, 6.7, 6.7, 6.7, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.7, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.8, 6.9, 6.9, 6.8, 6.8, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.7, 6.7, 6.7, 6.8, 6.8, 6.8, 6.9, 6.9, 7.0, 7.1, 7.1, 7.2, 7.2, 7.3, 7.3, 7.3, 7.4, 7.4, 7.4, 7.5, 7.6, 7.7, 7.8, 7.8, 7.9, 7.9, 8.0 },
    { 0.0, 0.0, 0.0, 6.0, 6.0, 6.2, 6.4, 6.4, 6.4, 6.4, 6.3, 6.2, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.5, 6.5, 6.5, 6.5, 6.6, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.6, 6.6, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.6, 6.6, 6.6, 6.6, 6.6, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.6, 6.6, 6.6, 6.7, 6.8, 6.8, 6.8, 6.9, 6.9, 7.0, 7.0, 7.1, 7.2, 7.2, 7.3, 7.3, 7.3, 7.4, 7.5, 7.6, 7.7, 7.8, 7.9, 7.9, 8.0, 8.0 },
    { 0.0, 0.0, 0.0, 5.9, 6.0, 6.0, 6.0, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.4, 6.4, 6.4, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.4, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.6, 6.6, 6.6, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.6, 6.6, 6.7, 6.7, 6.8, 6.8, 6.9, 6.9, 7.0, 7.1, 7.1, 7.2, 7.2, 7.2, 7.3, 7.4, 7.5, 7.7, 7.8, 7.9, 7.9, 8.0, 8.0 },
    { 0.0, 0.0, 0.0, 5.7, 5.7, 5.8, 5.8, 5.8, 5.9, 5.9, 6.0, 6.0, 6.1, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.6, 6.6, 6.6, 6.6, 6.7, 6.7, 6.7, 6.8, 6.8, 6.9, 7.0, 7.0, 7.1, 7.1, 7.2, 7.2, 7.3, 7.3, 7.4, 7.6, 7.8, 7.9, 8.0, 8.0, 8.0 },
    { 0.0, 0.0, 0.0, 5.6, 5.7, 5.7, 5.8, 5.8, 5.9, 6.0, 6.0, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.1, 6.1, 6.0, 6.0, 6.0, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.4, 6.4, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.8, 6.9, 6.9, 7.0, 7.0, 7.1, 7.2, 7.2, 7.3, 7.3, 7.4, 7.6, 7.8, 7.9, 8.0, 8.0, 8.0 },
    { 0.0, 0.0, 0.0, 5.7, 5.7, 5.8, 5.9, 6.0, 6.1, 6.2, 6.2, 6.2, 6.2, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.4, 6.4, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.7, 6.7, 6.7, 6.7, 6.7, 6.7, 6.8, 6.8, 6.9, 6.9, 7.0, 7.0, 7.1, 7.1, 7.2, 7.3, 7.3, 7.4, 7.6, 7.8, 7.9, 8.0, 8.0, 8.0 },
    { 0.0, 0.0, 0.0, 5.7, 5.8, 5.9, 6.0, 6.1, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.4, 6.5, 6.5, 6.5, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.6, 6.6, 6.6, 6.6, 6.6, 6.7, 6.7, 6.8, 6.9, 6.9, 6.9, 7.0, 7.1, 7.1, 7.2, 7.3, 7.4, 7.5, 7.6, 7.7, 7.8, 7.9, 8.0, 8.0 },
    { 0.0, 0.0, 0.0, 5.8, 5.9, 5.9, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.5, 6.5, 6.5, 6.5, 6.5, 6.6, 6.6, 6.7, 6.7, 6.8, 6.8, 6.9, 6.9, 7.0, 7.0, 7.1, 7.2, 7.2, 7.3, 7.4, 7.5, 7.6, 7.7, 7.8, 7.8, 7.8, 7.9 },
    { 0.0, 0.0, 0.0, 5.8, 5.8, 5.9, 6.0, 6.1, 6.2, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.2, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.4, 6.4, 6.5, 6.5, 6.6, 6.6, 6.7, 6.7, 6.8, 6.8, 6.9, 6.9, 7.0, 7.1, 7.1, 7.2, 7.2, 7.3, 7.4, 7.5, 7.6, 7.7, 7.7, 7.7, 7.8, 7.8 },
    { 0.0, 0.0, 0.0, 5.8, 5.8, 5.9, 5.9, 6.0, 6.1, 6.1, 6.2, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.4, 6.4, 6.5, 6.5, 6.6, 6.7, 6.7, 6.8, 6.8, 6.9, 6.9, 7.0, 7.0, 7.1, 7.1, 7.2, 7.2, 7.3, 7.4, 7.5, 7.6, 7.6, 7.7, 7.7, 7.7, 7.7 },
    { 0.0, 0.0, 0.0, 5.7, 5.8, 5.8, 5.8, 5.9, 5.8, 5.9, 6.0, 6.0, 6.1, 6.1, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.1, 6.1, 6.1, 6.0, 6.0, 6.0, 6.0, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.4, 6.4, 6.5, 6.5, 6.6, 6.7, 6.7, 6.8, 6.8, 6.8, 7.0, 7.0, 7.0, 7.1, 7.1, 7.2, 7.2, 7.3, 7.3, 7.4, 7.5, 7.6, 7.6, 7.6, 7.7, 7.7 },
    { 0.0, 0.0, 0.0, 5.7, 5.7, 5.7, 5.8, 5.8, 5.8, 5.8, 5.8, 5.9, 5.9, 5.9, 6.0, 6.1, 6.1, 6.1, 6.2, 6.2, 6.3, 6.3, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.4, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.1, 6.1, 6.1, 6.1, 6.1, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.4, 6.4, 6.4, 6.5, 6.5, 6.5, 6.6, 6.7, 6.7, 6.8, 6.9, 6.9, 7.0, 7.0, 7.0, 7.1, 7.1, 7.2, 7.2, 7.3, 7.3, 7.4, 7.5, 7.5, 7.5, 7.6, 7.6 },
    { 0.0, 0.0, 0.0, 5.7, 5.7, 5.7, 5.7, 5.7, 5.7, 5.8, 5.8, 5.8, 5.8, 5.8, 5.9, 5.9, 6.0, 6.0, 6.0, 6.0, 6.1, 6.1, 6.1, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.3, 6.4, 6.4, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.2, 6.2, 6.2, 6.2, 6.2, 6.2, 6.1, 6.1, 6.1, 6.1, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.1, 6.1, 6.1, 6.1, 6.2, 6.2, 6.2, 6.2, 6.2, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.3, 6.4, 6.4, 6.4, 6.4, 6.5, 6.5, 6.6, 6.7, 6.7, 6.8, 6.9, 6.9, 7.0, 7.0, 7.0, 7.1, 7.1, 7.2, 7.2, 7.3, 7.4, 7.4, 7.4, 7.5, 7.5 }
};

*/

/* Preliminary, unpublished Q function of Saul & Bormann (2007) */
static float __qmb[17][108] = {
    { // 0 km
0.0, 0.0, 0.0, 5.75, 5.88, 5.96, 6.01, 6.10, 6.14, 6.17, 6.17, 6.16, 6.15, 6.15, 6.10, 6.06, 6.05, 6.04, 6.05, 6.08, 6.11, 6.16, 6.23, 6.30, 6.39, 6.48, 6.55, 6.61, 6.64, 6.66, 6.67, 6.67, 6.67, 6.67, 6.66, 6.66, 6.65, 6.66, 6.66, 6.67, 6.67, 6.69, 6.69, 6.70, 6.71, 6.72, 6.72, 6.73, 6.74, 6.75, 6.76, 6.77, 6.77, 6.78, 6.79, 6.79, 6.80, 6.81, 6.82, 6.82, 6.83, 6.83, 6.84, 6.84, 6.84, 6.85, 6.85, 6.85, 6.86, 6.87, 6.87, 6.87, 6.87, 6.88, 6.88, 6.88, 6.89, 6.89, 6.90, 6.89, 6.89, 6.89, 6.90, 6.92, 6.94, 6.97, 7.01, 7.05, 7.09, 7.15, 7.20, 7.25, 7.31, 7.37, 7.42, 7.46, 7.50, 7.55, 7.59, 7.64, 7.68, 7.73, 7.75, 7.78, 7.78, 7.78, 7.78, 7.78
    },
    { // 25 km
0.0, 0.0, 0.0, 5.77, 5.85, 5.92, 5.99, 6.09, 6.16, 6.21, 6.23, 6.22, 6.20, 6.19, 6.15, 6.12, 6.10, 6.11, 6.12, 6.14, 6.17, 6.21, 6.26, 6.33, 6.41, 6.50, 6.57, 6.63, 6.65, 6.68, 6.69, 6.69, 6.69, 6.69, 6.69, 6.69, 6.69, 6.69, 6.68, 6.68, 6.68, 6.68, 6.69, 6.70, 6.71, 6.72, 6.74, 6.74, 6.74, 6.75, 6.75, 6.76, 6.77, 6.77, 6.78, 6.78, 6.78, 6.78, 6.79, 6.79, 6.80, 6.81, 6.81, 6.81, 6.82, 6.82, 6.83, 6.83, 6.83, 6.83, 6.83, 6.83, 6.83, 6.84, 6.84, 6.84, 6.85, 6.85, 6.85, 6.86, 6.86, 6.86, 6.86, 6.88, 6.90, 6.94, 6.98, 7.02, 7.06, 7.12, 7.16, 7.21, 7.26, 7.31, 7.36, 7.41, 7.45, 7.50, 7.54, 7.58, 7.62, 7.65, 7.67, 7.68, 7.68, 7.68, 7.68, 7.68
    },
    { // 50 km
0.0, 0.0, 0.0, 5.70, 5.86, 5.92, 5.96, 6.03, 6.10, 6.12, 6.15, 6.17, 6.16, 6.14, 6.13, 6.10, 6.08, 6.10, 6.12, 6.14, 6.16, 6.21, 6.24, 6.32, 6.40, 6.49, 6.54, 6.60, 6.63, 6.66, 6.69, 6.70, 6.69, 6.68, 6.68, 6.65, 6.65, 6.66, 6.66, 6.65, 6.66, 6.65, 6.65, 6.66, 6.68, 6.69, 6.71, 6.73, 6.74, 6.75, 6.76, 6.75, 6.74, 6.74, 6.75, 6.75, 6.76, 6.77, 6.77, 6.76, 6.76, 6.76, 6.76, 6.76, 6.78, 6.77, 6.79, 6.78, 6.78, 6.78, 6.79, 6.78, 6.78, 6.79, 6.79, 6.78, 6.78, 6.78, 6.78, 6.79, 6.79, 6.79, 6.80, 6.81, 6.83, 6.87, 6.92, 6.97, 7.02, 7.07, 7.12, 7.16, 7.21, 7.27, 7.33, 7.38, 7.41, 7.45, 7.49, 7.51, 7.54, 7.56, 7.58, 7.60, 7.60, 7.60, 7.60, 7.60
    },
    { // 75 km
0.0, 0.0, 0.0, 5.73, 5.80, 5.82, 5.88, 6.01, 6.10, 6.14, 6.15, 6.14, 6.14, 6.09, 6.04, 6.03, 6.03, 6.03, 6.07, 6.11, 6.13, 6.18, 6.22, 6.27, 6.37, 6.46, 6.53, 6.59, 6.62, 6.64, 6.66, 6.67, 6.67, 6.66, 6.65, 6.63, 6.64, 6.63, 6.63, 6.64, 6.63, 6.63, 6.64, 6.65, 6.66, 6.66, 6.67, 6.69, 6.72, 6.74, 6.75, 6.74, 6.74, 6.72, 6.72, 6.73, 6.72, 6.72, 6.73, 6.72, 6.73, 6.73, 6.74, 6.74, 6.74, 6.73, 6.74, 6.74, 6.74, 6.74, 6.75, 6.75, 6.76, 6.77, 6.78, 6.78, 6.77, 6.76, 6.76, 6.76, 6.74, 6.75, 6.77, 6.79, 6.82, 6.87, 6.93, 6.99, 7.03, 7.08, 7.13, 7.17, 7.21, 7.28, 7.34, 7.39, 7.41, 7.47, 7.53, 7.57, 7.62, 7.63, 7.64, 7.64, 7.64, 7.64, 7.64, 7.64
    },
    { // 100 km
0.0, 0.0, 0.0, 5.67, 5.73, 5.73, 5.83, 5.96, 6.02, 6.10, 6.20, 6.18, 6.15, 6.10, 6.05, 6.03, 6.01, 6.03, 6.05, 6.09, 6.10, 6.14, 6.19, 6.25, 6.33, 6.43, 6.50, 6.56, 6.59, 6.61, 6.61, 6.61, 6.60, 6.59, 6.57, 6.57, 6.57, 6.57, 6.56, 6.59, 6.58, 6.59, 6.59, 6.60, 6.59, 6.60, 6.60, 6.60, 6.61, 6.63, 6.65, 6.66, 6.67, 6.67, 6.67, 6.67, 6.68, 6.68, 6.68, 6.66, 6.66, 6.66, 6.66, 6.66, 6.66, 6.66, 6.65, 6.65, 6.66, 6.66, 6.68, 6.68, 6.69, 6.69, 6.69, 6.70, 6.70, 6.70, 6.71, 6.71, 6.71, 6.72, 6.73, 6.75, 6.77, 6.82, 6.86, 6.92, 6.98, 7.03, 7.07, 7.12, 7.17, 7.22, 7.27, 7.32, 7.34, 7.38, 7.42, 7.47, 7.51, 7.54, 7.54, 7.56, 7.56, 7.56, 7.56, 7.56
    }, 
    { // 150 km
0.0, 0.0, 0.0, 5.64, 5.68, 5.70, 5.78, 5.90, 5.98, 6.03, 6.08, 6.03, 5.94, 5.88, 5.91, 5.92, 5.95, 6.00, 6.03, 6.04, 6.06, 6.09, 6.12, 6.18, 6.24, 6.33, 6.41, 6.47, 6.52, 6.51, 6.51, 6.49, 6.48, 6.45, 6.45, 6.45, 6.46, 6.46, 6.46, 6.47, 6.46, 6.46, 6.46, 6.47, 6.47, 6.50, 6.50, 6.52, 6.51, 6.52, 6.53, 6.53, 6.52, 6.54, 6.55, 6.55, 6.56, 6.56, 6.56, 6.55, 6.56, 6.56, 6.57, 6.56, 6.57, 6.57, 6.56, 6.56, 6.57, 6.57, 6.58, 6.59, 6.59, 6.59, 6.60, 6.60, 6.61, 6.62, 6.63, 6.63, 6.64, 6.65, 6.66, 6.68, 6.71, 6.75, 6.78, 6.84, 6.89, 6.94, 6.99, 7.04, 7.10, 7.14, 7.20, 7.27, 7.30, 7.33, 7.38, 7.42, 7.46, 7.49, 7.51, 7.53, 7.53, 7.53, 7.53, 7.53
    },
    { // 200 km
0.0, 0.0, 0.0, 5.49, 5.50, 5.56, 5.68, 5.76, 5.86, 5.85, 5.84, 5.73, 5.70, 5.68, 5.75, 5.77, 5.84, 5.89, 5.96, 5.94, 5.97, 6.01, 6.07, 6.13, 6.22, 6.29, 6.34, 6.38, 6.42, 6.41, 6.41, 6.40, 6.38, 6.35, 6.35, 6.36, 6.36, 6.37, 6.38, 6.39, 6.38, 6.37, 6.35, 6.34, 6.35, 6.37, 6.38, 6.41, 6.40, 6.41, 6.42, 6.42, 6.42, 6.44, 6.44, 6.42, 6.44, 6.44, 6.44, 6.46, 6.48, 6.47, 6.47, 6.45, 6.43, 6.43, 6.42, 6.43, 6.46, 6.47, 6.47, 6.49, 6.50, 6.51, 6.52, 6.53, 6.53, 6.54, 6.55, 6.54, 6.55, 6.55, 6.56, 6.58, 6.61, 6.64, 6.69, 6.73, 6.76, 6.81, 6.85, 6.89, 6.95, 7.03, 7.08, 7.17, 7.21, 7.22, 7.32, 7.35, 7.42, 7.45, 7.47, 7.47, 7.47, 7.47, 7.47, 7.47
    },
    { // 250 km
0.0, 0.0, 0.0, 5.49, 5.49, 5.57, 5.73, 5.82, 5.90, 5.85, 5.76, 5.66, 5.66, 5.64, 5.73, 5.76, 5.87, 5.88, 5.95, 5.94, 5.98, 6.02, 6.09, 6.15, 6.22, 6.26, 6.30, 6.35, 6.38, 6.38, 6.37, 6.35, 6.31, 6.28, 6.28, 6.29, 6.28, 6.30, 6.30, 6.28, 6.28, 6.28, 6.28, 6.30, 6.32, 6.32, 6.33, 6.34, 6.31, 6.33, 6.33, 6.35, 6.34, 6.38, 6.37, 6.38, 6.39, 6.38, 6.38, 6.40, 6.43, 6.43, 6.46, 6.43, 6.41, 6.39, 6.39, 6.39, 6.41, 6.44, 6.44, 6.45, 6.46, 6.47, 6.47, 6.50, 6.51, 6.53, 6.54, 6.56, 6.56, 6.56, 6.57, 6.59, 6.60, 6.63, 6.68, 6.73, 6.75, 6.80, 6.84, 6.89, 6.95, 7.03, 7.09, 7.17, 7.24, 7.27, 7.31, 7.35, 7.42, 7.43, 7.44, 7.47, 7.47, 7.47, 7.47, 7.47
    },
    { // 300 km
0.0, 0.0, 0.0, 5.52, 5.47, 5.55, 5.68, 5.76, 5.82, 5.84, 5.77, 5.70, 5.72, 5.72, 5.77, 5.77, 5.83, 5.84, 5.87, 5.88, 5.93, 5.99, 6.09, 6.20, 6.24, 6.27, 6.30, 6.31, 6.32, 6.32, 6.30, 6.28, 6.24, 6.23, 6.21, 6.21, 6.21, 6.23, 6.22, 6.22, 6.23, 6.22, 6.23, 6.24, 6.24, 6.25, 6.25, 6.23, 6.21, 6.24, 6.25, 6.28, 6.30, 6.34, 6.31, 6.34, 6.33, 6.31, 6.32, 6.35, 6.36, 6.38, 6.42, 6.41, 6.42, 6.42, 6.40, 6.39, 6.41, 6.43, 6.41, 6.42, 6.44, 6.44, 6.45, 6.48, 6.50, 6.52, 6.53, 6.54, 6.53, 6.53, 6.52, 6.55, 6.57, 6.61, 6.66, 6.70, 6.72, 6.76, 6.80, 6.84, 6.90, 6.97, 7.03, 7.11, 7.18, 7.22, 7.27, 7.30, 7.36, 7.38, 7.41, 7.42, 7.42, 7.42, 7.42, 7.42
    },
    { // 350 km
0.0, 0.0, 0.0, 5.58, 5.51, 5.55, 5.64, 5.74, 5.75, 5.77, 5.74, 5.72, 5.75, 5.78, 5.79, 5.79, 5.82, 5.82, 5.83, 5.87, 5.93, 5.99, 6.06, 6.15, 6.16, 6.18, 6.23, 6.27, 6.25, 6.27, 6.25, 6.21, 6.19, 6.21, 6.20, 6.20, 6.19, 6.18, 6.16, 6.17, 6.15, 6.16, 6.19, 6.20, 6.18, 6.21, 6.20, 6.20, 6.17, 6.20, 6.20, 6.24, 6.25, 6.31, 6.29, 6.32, 6.31, 6.29, 6.29, 6.30, 6.30, 6.35, 6.38, 6.38, 6.42, 6.42, 6.39, 6.40, 6.40, 6.40, 6.38, 6.37, 6.38, 6.38, 6.40, 6.43, 6.46, 6.49, 6.50, 6.52, 6.51, 6.51, 6.52, 6.56, 6.58, 6.61, 6.66, 6.71, 6.73, 6.77, 6.81, 6.85, 6.91, 6.97, 7.03, 7.09, 7.17, 7.22, 7.27, 7.30, 7.35, 7.36, 7.38, 7.38, 7.38, 7.38, 7.38, 7.38
    },
    { // 400 km
0.0, 0.0, 0.0, 5.44, 5.44, 5.46, 5.51, 5.61, 5.61, 5.71, 5.72, 5.74, 5.76, 5.79, 5.78, 5.78, 5.78, 5.80, 5.81, 5.86, 5.88, 5.99, 6.05, 6.15, 6.15, 6.16, 6.18, 6.21, 6.18, 6.19, 6.19, 6.17, 6.14, 6.17, 6.15, 6.14, 6.14, 6.14, 6.12, 6.14, 6.13, 6.13, 6.16, 6.18, 6.15, 6.16, 6.16, 6.17, 6.14, 6.18, 6.18, 6.21, 6.22, 6.28, 6.27, 6.28, 6.28, 6.26, 6.27, 6.26, 6.27, 6.30, 6.32, 6.32, 6.38, 6.38, 6.36, 6.38, 6.38, 6.37, 6.36, 6.37, 6.37, 6.39, 6.40, 6.42, 6.44, 6.46, 6.46, 6.47, 6.48, 6.49, 6.49, 6.53, 6.57, 6.62, 6.66, 6.71, 6.74, 6.77, 6.81, 6.86, 6.92, 6.99, 7.06, 7.12, 7.20, 7.24, 7.29, 7.32, 7.36, 7.37, 7.41, 7.39, 7.39, 7.39, 7.39, 7.39
    },
    { // 450 km
0.0, 0.0, 0.0, 5.38, 5.43, 5.49, 5.52, 5.57, 5.59, 5.66, 5.73, 5.76, 5.77, 5.78, 5.79, 5.81, 5.82, 5.84, 5.85, 5.89, 5.95, 6.02, 6.06, 6.11, 6.15, 6.18, 6.18, 6.17, 6.17, 6.18, 6.16, 6.16, 6.15, 6.14, 6.15, 6.14, 6.14, 6.14, 6.14, 6.14, 6.14, 6.16, 6.15, 6.14, 6.14, 6.14, 6.14, 6.16, 6.17, 6.18, 6.19, 6.21, 6.23, 6.23, 6.24, 6.25, 6.26, 6.26, 6.26, 6.27, 6.27, 6.29, 6.31, 6.31, 6.32, 6.33, 6.34, 6.35, 6.34, 6.35, 6.35, 6.36, 6.37, 6.38, 6.39, 6.41, 6.41, 6.42, 6.44, 6.44, 6.46, 6.47, 6.49, 6.52, 6.54, 6.58, 6.62, 6.66, 6.69, 6.73, 6.78, 6.81, 6.87, 6.92, 6.99, 7.04, 7.10, 7.14, 7.19, 7.22, 7.25, 7.27, 7.28, 7.32, 7.32, 7.32, 7.32, 7.32
    },
    { // 500 km
0.0, 0.0, 0.0, 5.37, 5.40, 5.44, 5.48, 5.54, 5.57, 5.65, 5.71, 5.75, 5.75, 5.78, 5.78, 5.80, 5.82, 5.85, 5.85, 5.90, 5.94, 6.00, 6.04, 6.08, 6.12, 6.14, 6.14, 6.14, 6.14, 6.15, 6.15, 6.14, 6.13, 6.13, 6.14, 6.14, 6.12, 6.11, 6.11, 6.11, 6.11, 6.12, 6.11, 6.11, 6.11, 6.12, 6.13, 6.15, 6.16, 6.17, 6.18, 6.20, 6.21, 6.21, 6.22, 6.23, 6.23, 6.23, 6.23, 6.25, 6.25, 6.26, 6.30, 6.29, 6.31, 6.32, 6.34, 6.34, 6.34, 6.34, 6.34, 6.35, 6.36, 6.37, 6.37, 6.38, 6.38, 6.39, 6.41, 6.42, 6.44, 6.45, 6.47, 6.50, 6.52, 6.57, 6.61, 6.65, 6.68, 6.72, 6.77, 6.81, 6.86, 6.93, 6.99, 7.04, 7.10, 7.13, 7.18, 7.20, 7.23, 7.25, 7.28, 7.34, 7.34, 7.34, 7.34, 7.34
    },
    { // 550 km
0.0, 0.0, 0.0, 5.42, 5.44, 5.48, 5.50, 5.55, 5.57, 5.63, 5.66, 5.70, 5.72, 5.75, 5.78, 5.80, 5.84, 5.88, 5.91, 5.95, 5.98, 6.05, 6.10, 6.13, 6.15, 6.16, 6.15, 6.14, 6.12, 6.13, 6.12, 6.12, 6.11, 6.11, 6.11, 6.10, 6.10, 6.10, 6.10, 6.10, 6.09, 6.10, 6.09, 6.09, 6.10, 6.11, 6.13, 6.15, 6.17, 6.19, 6.20, 6.22, 6.22, 6.21, 6.21, 6.21, 6.22, 6.21, 6.21, 6.23, 6.23, 6.25, 6.28, 6.28, 6.29, 6.30, 6.31, 6.30, 6.30, 6.31, 6.30, 6.31, 6.32, 6.33, 6.34, 6.35, 6.37, 6.38, 6.40, 6.41, 6.43, 6.45, 6.47, 6.50, 6.53, 6.57, 6.61, 6.65, 6.68, 6.73, 6.77, 6.81, 6.86, 6.92, 6.96, 7.01, 7.07, 7.11, 7.15, 7.18, 7.20, 7.23, 7.25, 7.29, 7.29, 7.29, 7.29, 7.29
    },
    { // 600 km
0.0, 0.0, 0.0, 5.43, 5.44, 5.48, 5.50, 5.53, 5.53, 5.59, 5.61, 5.65, 5.69, 5.74, 5.78, 5.82, 5.88, 5.94, 6.01, 6.02, 6.08, 6.15, 6.20, 6.22, 6.21, 6.21, 6.20, 6.16, 6.13, 6.13, 6.11, 6.11, 6.12, 6.12, 6.11, 6.10, 6.10, 6.11, 6.11, 6.10, 6.11, 6.09, 6.10, 6.10, 6.10, 6.12, 6.14, 6.16, 6.17, 6.18, 6.19, 6.21, 6.21, 6.20, 6.19, 6.19, 6.19, 6.19, 6.20, 6.21, 6.21, 6.22, 6.25, 6.26, 6.27, 6.27, 6.28, 6.26, 6.27, 6.27, 6.26, 6.26, 6.28, 6.29, 6.30, 6.31, 6.34, 6.36, 6.37, 6.39, 6.42, 6.45, 6.48, 6.51, 6.53, 6.56, 6.60, 6.64, 6.68, 6.72, 6.76, 6.80, 6.85, 6.91, 6.96, 7.01, 7.07, 7.11, 7.14, 7.18, 7.19, 7.22, 7.25, 7.28, 7.28, 7.28, 7.28, 7.28
    },
    { // 650 km
0.0, 0.0, 0.0, 5.41, 5.31, 5.38, 5.37, 5.37, 5.43, 5.45, 5.51, 5.58, 5.64, 5.69, 5.76, 5.85, 5.89, 5.98, 6.10, 6.09, 6.19, 6.24, 6.23, 6.21, 6.20, 6.20, 6.16, 6.11, 6.10, 6.07, 6.03, 6.05, 6.10, 6.10, 6.10, 6.13, 6.14, 6.15, 6.13, 6.14, 6.15, 6.14, 6.15, 6.14, 6.13, 6.15, 6.16, 6.18, 6.17, 6.18, 6.18, 6.19, 6.20, 6.18, 6.19, 6.19, 6.19, 6.20, 6.20, 6.21, 6.21, 6.21, 6.25, 6.24, 6.26, 6.27, 6.28, 6.28, 6.28, 6.27, 6.26, 6.25, 6.27, 6.28, 6.29, 6.29, 6.30, 6.32, 6.33, 6.35, 6.39, 6.43, 6.47, 6.49, 6.52, 6.55, 6.58, 6.64, 6.69, 6.73, 6.77, 6.82, 6.88, 6.92, 6.96, 7.03, 7.08, 7.10, 7.12, 7.16, 7.20, 7.21, 7.25, 7.28, 7.28, 7.28, 7.28, 7.28
    },
    { // 700 km
0.0, 0.0, 0.0, 4.89, 5.10, 5.29, 5.29, 5.27, 5.27, 5.44, 5.64, 5.67, 5.67, 5.79, 5.79, 5.89, 6.17, 6.05, 6.14, 6.16, 6.13, 6.16, 6.11, 6.08, 6.09, 6.07, 6.04, 6.00, 6.03, 6.01, 6.03, 6.05, 6.05, 6.08, 6.12, 6.15, 6.17, 6.21, 6.23, 6.21, 6.22, 6.19, 6.17, 6.17, 6.16, 6.16, 6.12, 6.10, 6.09, 6.10, 6.12, 6.11, 6.10, 6.11, 6.11, 6.11, 6.09, 6.10, 6.19, 6.20, 6.21, 6.20, 6.23, 6.20, 6.21, 6.19, 6.23, 6.22, 6.22, 6.20, 6.20, 6.23, 6.24, 6.27, 6.28, 6.33, 6.35, 6.39, 6.42, 6.47, 6.53, 6.54, 6.58, 6.58, 6.60, 6.62, 6.64, 6.69, 6.70, 6.75, 6.74, 6.76, 6.83, 6.88, 6.96, 6.96, 7.00, 7.03, 7.07, 7.13, 7.18, 7.18, 7.21, 7.24, 7.24, 7.24, 7.24, 7.24
    }
};

#define max(x,y) ((x)>(y)?(x):(y))
#define min(x,y) ((x)<(y)?(x):(y))

static double
bmagnz(double amplitude, double period, double delta, double depth, int *err)
{
	int j, k;
	double mb, s1, s2, q1, Q;

	*err = 0;
	
	if (depth < 0. || depth > 700. || delta < 5 || delta > 105) {
		*err = -1;
		return 0.;
	}
	
	if (depth < 100) {
		k  = (int) (depth/25.) + 2;
		s1 = 0.04*(depth-25.*(k-2));
	}
	else {
		k = min((int)(depth/50.)+4, 17);
		s1 = .02*(depth-50.*(k-4));
	}
	j  = max(min((int)(delta), 108), 2);
	s2 = delta-j;
	q1 = __qmb[k-2][j-2] + s1*(__qmb[k-1][j-2]-__qmb[k-2][j-2]);
	Q  = q1+s2*((__qmb[k-2][j-1]+s1*(__qmb[k-1][j-1]-__qmb[k-2][j-1]))-q1);
	
	mb = log10(amplitude/period) + Q;

	return mb;
}

namespace Seiscomp {
namespace Magnitudes {

bool
compute_mb(
	double amplitude, // in micrometers
	double period,    // in seconds
	double delta,     // in degrees
	double depth,     // in kilometers
	double *mag)
{
	int err;
	double m = bmagnz(amplitude, period, delta, depth, &err);
	if (err) return false;

	*mag = m;

	return true;
}

bool
compute_mb_fromVelocity(
	double amplitude, // in micrometers/second
	double delta,     // in degrees
	double depth,     // in kilometers
	double *mag)
{
	int err;
	double m = bmagnz(amplitude, 2*M_PI, delta, depth, &err);
	if (err) return false;

	*mag = m;

	return true;
}

} // namespace Magnitudes
} // namespace Seiscomp
