"""Paired-read decision rules with exact binomial tails and abstention."""
from dataclasses import dataclass
from fractions import Fraction as F


def binomial_weights(n,p):
    p=F(p)
    if not isinstance(n,int) or n<0 or not 0<=p<=1:raise ValueError('Integer sample size and probability required.')
    a,d=p.numerator,p.denominator;b=d-a
    if p==1:return [0]*n+[1],1
    if p==0:return [1]+[0]*n,1
    w=[b**n]
    for j in range(n):
        num=w[-1]*(n-j)*a;den=(j+1)*b;q,r=divmod(num,den)
        if r:raise ArithmeticError('Nonintegral binomial recurrence.')
        w.append(q)
    if sum(w)!=d**n:raise ArithmeticError('Binomial mass mismatch.')
    return w,d**n


def tail(n,p,count,upper=False):
    if not isinstance(count,int) or not 0<=count<=n:raise ValueError('Count outside sample.')
    w,d=binomial_weights(n,p);return F(sum(w[count:] if upper else w[:count+1]),d)


@dataclass(frozen=True)
class SymmetricReadout:
    flip_probability:F=F(1,10)
    def __post_init__(self):
        object.__setattr__(self,'flip_probability',F(self.flip_probability))
        if not 0<=self.flip_probability<=1:raise ValueError('Flip probability outside [0,1].')
    @property
    def attenuation(self):return (1-2*self.flip_probability)**2
    def agreement(self,c):
        c=F(c)
        if not -F(1,4)<=c<=F(1,4):raise ValueError('Covariance outside source class.')
        return F(1,2)+2*self.attenuation*c


class PairedDecisionRule:
    def __init__(self,band,samples,readout=SymmetricReadout(),alpha=F(1,20)):
        self.band=band;self.samples=samples;self.readout=readout;self.alpha=F(alpha)
        if not isinstance(samples,int) or samples<1 or not 0<self.alpha<1:raise ValueError('Positive sample size and alpha in (0,1) required.')
        if readout.attenuation==0:raise ValueError('A pure-noise channel supplies no calibrated covariance information.')
        self.low_cut=-1;self.high_cut=samples+1
        self.p_minus=F(1,2)+2*readout.attenuation*band.low_boundary
        self.p_plus=F(1,2)+2*readout.attenuation*band.high_boundary
        if band.low_boundary>F(1,4):self.low_cut=samples
        elif band.low_boundary>=-F(1,4):
            w,d=binomial_weights(samples,self.p_minus);cum=0
            for j,a in enumerate(w):
                cum+=a
                if 2*cum*self.alpha.denominator<=self.alpha.numerator*d:self.low_cut=j
        if band.high_boundary<-F(1,4):self.high_cut=0
        elif band.high_boundary<=F(1,4):
            w,d=binomial_weights(samples,self.p_plus);remaining=d
            for j,a in enumerate(w):
                if 2*remaining*self.alpha.denominator<=self.alpha.numerator*d:self.high_cut=j;break
                remaining-=a
        if self.low_cut>=self.high_cut:raise ArithmeticError('Decision regions overlap.')
    def decide(self,count):
        if not isinstance(count,int) or not 0<=count<=self.samples:raise ValueError('Agreement count outside sample.')
        return self.band.low_action if count<=self.low_cut else self.band.high_action if count>=self.high_cut else 'unresolved'
    def probabilities(self,c):
        w,d=binomial_weights(self.samples,self.readout.agreement(c));low=F(sum(w[:self.low_cut+1]),d);high=F(sum(w[self.high_cut:]),d)
        return dict(low=low,high=high,unresolved=1-low-high,resolution=low+high)


def binomial_confidence_interval(n,count,alpha=F(1,20),bits=32):
    """Outward rational bisection of equal-tailed Clopper-Pearson endpoints."""
    alpha=F(alpha)
    if not isinstance(n,int) or n<1 or not isinstance(count,int) or not 0<=count<=n or not 0<alpha<1 or bits<1:raise ValueError('Invalid binomial interval request.')
    result=[]
    for is_lower in [True,False]:
        if is_lower and count==0:result.append(F(0));continue
        if not is_lower and count==n:result.append(F(1));continue
        left,right=F(0),F(1)
        for _ in range(bits):
            mid=(left+right)/2;prob=tail(n,mid,count,upper=is_lower)
            if (prob<alpha/2)==is_lower:left=mid
            else:right=mid
        result.append(left if is_lower else right)
    return tuple(result)


def project_covariance(bit_interval,attenuation_interval):
    lo,hi=map(F,bit_interval);kl,ku=map(F,attenuation_interval)
    if not 0<=lo<=hi<=1 or not 0<=kl<=ku<=1:raise ValueError('Invalid bit or calibration interval.')
    if kl==0:
        if lo<=F(1,2)<=hi:return (-F(1,4),F(1,4))
        raise ValueError('Calibration includes zero attenuation without a compatible centered bit value; refine the joint feasible-set calculation.')
    candidates=[(p-F(1,2))/(2*k) for p in [lo,hi] for k in [kl,ku]]
    low,high=max(-F(1,4),min(candidates)),min(F(1,4),max(candidates))
    if low>high:raise ValueError('Empty joint source/calibration set.')
    return low,high


def calibrated_free_pulse_choice(band,bit_interval):
    lo,hi=map(F,bit_interval)
    if band.low_action!='AB' or band.high_action!='BA' or band.low_boundary<=0:raise ValueError('Calibration-free rule requires the positive pulse margin.')
    if not 0<=lo<=hi<=1:raise ValueError('Invalid confidence interval.')
    if hi<F(1,2):return 'AB'
    if lo>F(1,2)+2*band.high_boundary:return 'BA'
    return 'unresolved'


def covariance_sign(observed_interval,error_covariance_bound):
    lo,hi=map(F,observed_interval);delta=F(error_covariance_bound)
    if lo>hi or delta<0:raise ValueError('Ordered interval and nonnegative error reserve required.')
    return 'positive' if lo>delta else 'negative' if hi<-delta else 'unresolved'


def minimax_regret(risk_difference_interval):
    # R_first - R_second, explicitly NOT the extinction difference.
    L,U=map(F,risk_difference_interval)
    if L>U:raise ValueError('Ordered risk interval required.')
    if U<=0:return dict(probability_first=F(1),regret=F(0))
    if L>=0:return dict(probability_first=F(0),regret=F(0))
    return dict(probability_first=-L/(U-L),regret=-L*U/(U-L))