#ifndef CPC_CONFIDENCE_HPP_
#define CPC_CONFIDENCE_HPP_
#include <cmath>
#include <stdexcept>
#include "cpc_sketch.hpp"
namespace datasketches {
static const double ICON_ERROR_CONSTANT = 0.693147180559945286;
static const int16_t ICON_LOW_SIDE_DATA [33] = { 6037, 5720, 5328, 6411, 6262, 5682, 6724, 6403, 6127, 6665, 6411, 6208, 6959, 6525, 6427, 6892, 6665, 6619, 6792, 6752, 6690, 6899, 6818, 6708, 6871, 6845, 6812, 6909, 6861, 6828, 6919, 6897, 6842, };
static const int16_t ICON_HIGH_SIDE_DATA [33] = { 8031, 8559, 9309, 7084, 7959, 8660, 7141, 7514, 7876, 7458, 7430, 7572, 6892, 7141, 7497, 6889, 7132, 7290, 7075, 7118, 7185, 7040, 7047, 7085, 6993, 7019, 7053, 6953, 7001, 6983, 6944, 6966, 7004, };
static const double HIP_ERROR_CONSTANT = 0.588705011257737332;
static const int16_t HIP_LOW_SIDE_DATA [33] = { 5871, 5247, 4826, 5877, 5403, 5070, 5873, 5533, 5304, 5878, 5632, 5464, 5874, 5690, 5564, 5880, 5745, 5619, 5875, 5784, 5701, 5866, 5789, 5742, 5869, 5827, 5784, 5876, 5860, 5827, 5881, 5853, 5842, };
static const int16_t HIP_HIGH_SIDE_DATA [33] = { 5855, 6688, 7391, 5886, 6444, 6923, 5885, 6254, 6594, 5889, 6134, 6326, 5900, 6072, 6203, 5875, 6005, 6089, 5871, 5980, 6040, 5889, 5941, 6015, 5871, 5926, 5973, 5866, 5901, 5915, 5880, 5914, 5953, };
template<typename A>
double get_icon_confidence_lb(const cpc_sketch_alloc<A>& sketch, int kappa) {
if (sketch.get_num_coupons() == 0) return 0.0;
const int lg_k = sketch.get_lg_k();
const long k = 1 << lg_k;
if (lg_k < 4) throw std::logic_error("lgk < 4");
if (kappa < 1 || kappa > 3) throw std::invalid_argument("kappa must be between 1 and 3");
double x = ICON_ERROR_CONSTANT;
if (lg_k <= 14) x = ((double) ICON_HIGH_SIDE_DATA[3 * (lg_k - 4) + (kappa - 1)]) / 10000.0;
const double rel = x / sqrt(k);
const double eps = kappa * rel;
const double est = sketch.get_icon_estimate();
double result = est / (1.0 + eps);
const double check = sketch.get_num_coupons();
if (result < check) result = check;
return result;
}
template<typename A>
double get_icon_confidence_ub(const cpc_sketch_alloc<A>& sketch, int kappa) {
if (sketch.get_num_coupons() == 0) return 0.0;
const int lg_k = sketch.get_lg_k();
const long k = 1 << lg_k;
if (lg_k < 4) throw std::logic_error("lgk < 4");
if (kappa < 1 || kappa > 3) throw std::invalid_argument("kappa must be between 1 and 3");
double x = ICON_ERROR_CONSTANT;
if (lg_k <= 14) x = ((double) ICON_LOW_SIDE_DATA[3 * (lg_k - 4) + (kappa - 1)]) / 10000.0;
const double rel = x / sqrt(k);
const double eps = kappa * rel;
const double est = sketch.get_icon_estimate();
const double result = est / (1.0 - eps);
return ceil(result); }
template<typename A>
double get_hip_confidence_lb(const cpc_sketch_alloc<A>& sketch, int kappa) {
if (sketch.get_num_coupons() == 0) return 0.0;
const int lg_k = sketch.get_lg_k();
const long k = 1 << lg_k;
if (lg_k < 4) throw std::logic_error("lgk < 4");
if (kappa < 1 || kappa > 3) throw std::invalid_argument("kappa must be between 1 and 3");
double x = HIP_ERROR_CONSTANT;
if (lg_k <= 14) x = ((double) HIP_HIGH_SIDE_DATA[3 * (lg_k - 4) + (kappa - 1)]) / 10000.0;
const double rel = x / (sqrt((double) k));
const double eps = ((double) kappa) * rel;
const double est = sketch.get_hip_estimate();
double result = est / (1.0 + eps);
const double check = (double) sketch.get_num_coupons();
if (result < check) result = check;
return result;
}
template<typename A>
double get_hip_confidence_ub(const cpc_sketch_alloc<A>& sketch, int kappa) {
if (sketch.get_num_coupons() == 0) return 0.0;
const int lg_k = sketch.get_lg_k();
const long k = 1 << lg_k;
if (lg_k < 4) throw std::logic_error("lgk < 4");
if (kappa < 1 || kappa > 3) throw std::invalid_argument("kappa must be between 1 and 3");
double x = HIP_ERROR_CONSTANT;
if (lg_k <= 14) x = ((double) HIP_LOW_SIDE_DATA[3 * (lg_k - 4) + (kappa - 1)]) / 10000.0;
const double rel = x / sqrt(k);
const double eps = kappa * rel;
const double est = sketch.get_hip_estimate();
const double result = est / (1.0 - eps);
return ceil(result); }
}
#endif