"""Finite original counts, one latent state per specimen, independent references."""
from dataclasses import dataclass
from fractions import Fraction as F


def probability(x):
    x=F(x)
    if not 0<=x<=1: raise ValueError('Probability outside [0,1].')
    return x


def positive_integer(n):
    if type(n)!=int or n<1: raise ValueError('Positive integer required.')
    return n


@dataclass(frozen=True)
class SpecimenBudget:
    original_ml: F=F(1)
    available_ml: F=F(1)
    loss_per_preparation_ml: F=F(1,20)

    def __post_init__(self):
        if not 0<self.available_ml<=self.original_ml or self.loss_per_preparation_ml<0:
            raise ValueError('Require 0 < available <= original and nonnegative loss.')

    def fractions(self, split=F(1,2)):
        split=probability(split)
        a=(self.available_ml*split-self.loss_per_preparation_ml)/self.original_ml
        b=(self.available_ml*(1-split)-self.loss_per_preparation_ml)/self.original_ml
        if min(a,b)<0: raise ValueError('Preparation loss exceeds allocated input.')
        return a,b

    def single(self):
        s=(self.available_ml-self.loss_per_preparation_ml)/self.original_ml
        return probability(s)


@dataclass(frozen=True)
class SourceLaw:
    # Each atom is (probability, recovery in A, recovery in B).
    atoms: tuple

    def __post_init__(self):
        if not self.atoms or sum(F(t[0]) for t in self.atoms)!=1:
            raise ValueError('Source weights must sum exactly to one.')
        for atom in self.atoms:
            if len(atom)!=3: raise ValueError('Expected (weight,X,Y).')
            for x in atom: probability(x)

    @classmethod
    def corners(cls,p1,p2,q):
        p1,p2,q=map(F,(p1,p2,q))
        return cls(((1-p1-p2+q,F(0),F(0)),(p1-q,F(1),F(0)),
                    (p2-q,F(0),F(1)),(q,F(1),F(1))))

    def moments(self):
        return tuple(sum(F(t[0])*v(t) for t in self.atoms) for v in
                     (lambda t:F(t[1]),lambda t:F(t[2]),lambda t:F(t[1])*F(t[2])))

    def pair_outcomes(self):
        p1,p2,q=self.moments()
        return {'00':1-p1-p2+q,'10':p1-q,'01':p2-q,'11':q}

    def negative(self,n,a,b):
        if type(n)!=int or n<0: raise ValueError('Original count must be a nonnegative integer.')
        a,b=probability(a),probability(b)
        if a+b>1: raise ValueError('Aliquots cannot duplicate the specimen.')
        return sum(F(w)*(1-a*F(x)-b*F(y))**n for w,x,y in self.atoms)

    def joint_error(self,n,a,b,m,gate='at-least-one'):
        positive_integer(m);p=self.pair_outcomes()
        if gate not in ('at-least-one','exactly-one'): raise ValueError('Unknown gate.')
        w=p['10']+p['01']+(p['11'] if gate=='at-least-one' else 0)
        return w**m*self.negative(n,a,b)

    def sample_specimen(self,n,a,b,rng):
        # Draw state ONCE; multinomial allocates original targets without replacement
        # between mutually exclusive routes. Counts are not two Poisson aliquots.
        self.negative(n,a,b)
        k=rng.choice(len(self.atoms),p=[float(t[0]) for t in self.atoms])
        _,x,y=self.atoms[k];a,b=F(a),F(b);x,y=F(x),F(y)
        masses=[1-a-b,a*(1-x),a*x,b*(1-y),b*y]
        counts=rng.multinomial(n,[float(v) for v in masses])
        return dict(state=k,path_counts=counts.tolist(),detected=int(counts[2]+counts[4]))


def gate_bound(H,m):
    H=probability(H);positive_integer(m)
    z=F(1) if H==1 else min(F(1),F(m,m+1)/(1-H))
    return z**m*(1-(1-H)*z),z


def minimum_pairs(H,alpha,cap=10000):
    H,alpha=probability(H),probability(alpha);positive_integer(cap)
    if H>alpha or alpha==0: return {'status':'infeasible','pairs':None}
    for m in range(1,cap+1):
        v,_=gate_bound(H,m)
        if v<=alpha: return {'status':'certified','pairs':m,'bound':v}
    return {'status':'budget-exhausted','pairs':None}


def moment_bound(n,a,b,p1,p2,q):
    return SourceLaw.corners(p1,p2,q).negative(positive_integer(n),a,b)


@dataclass(frozen=True)
class ExclusionPolicy:
    a:F
    b:F
    count_threshold:int
    reference_pairs:int
    alpha:F=F(1,20)

    def __post_init__(self):
        probability(self.a);probability(self.b);probability(self.alpha)
        if self.a+self.b>1: raise ValueError('Aliquot fractions exceed original specimen.')
        positive_integer(self.count_threshold);positive_integer(self.reference_pairs)

    def bound(self):
        # Unrestricted means: adversary uses the smaller aliquot.
        return gate_bound(max((1-self.a)**self.count_threshold,(1-self.b)**self.count_threshold),self.reference_pairs)[0]

    def report(self,reference_outcomes,target_detected):
        # Missing/invalid outcomes cannot silently disappear from the denominator.
        if len(reference_outcomes)!=self.reference_pairs or any(v not in ('00','10','01','11') for v in reference_outcomes):
            return {'issued':False,'reason':'incomplete or invalid reference record'}
        if type(target_detected)!=int or target_detected<0:
            return {'issued':False,'reason':'invalid target record'}
        if self.bound()>self.alpha: return {'issued':False,'reason':'design bound exceeds allowance'}
        if '00' in reference_outcomes or target_detected!=0:
            return {'issued':False,'reason':'gate failed or target detected'}
        return {'issued':True,'statement':f'fewer than {self.count_threshold} original targets',
                'joint_false_exclusion_upper':self.bound(),
                'scope':'Prespecified single experiment under the model; not a posterior probability or proof of zero targets.'}