"""Reproduce synthetic HPLC injection sequences. Python3.9+, standard library only."""
from decimal import Decimal, localcontext
from pathlib import Path
import argparse
import csv
import math

D = Decimal
RULES = {
    'resolution_min': '1.50',
    'blank_window_start_min': '3.80',
    'blank_window_end_min': '4.20',
    'blank_apex_flag_au': '0.020',
    'trace_start_min': '0.00',
    'trace_end_min': '8.00',
    'trace_step_min': '0.02',
}
# An empty list means no modeled components, not a detection-limit result.
CASES = {
    'A': {'blank': [], 'reference': [('R1','4.00','0.20','1.00'),('R2','5.60','0.30','0.60')],
          'sample': [('S1','4.00','0.20','0.80'),('S2','5.60','0.30','0.30')]},
    'B': {'blank': [], 'reference': [('R1','4.00','0.20','1.00'),('R2','5.20','0.30','0.60')],
          'sample': [('S1','4.00','0.20','0.80'),('S2','5.20','0.30','0.30')]},
    'C': {'blank': [('U1','4.02','0.20','0.060')],
          'reference': [('R1','4.00','0.20','1.00'),('R2','5.60','0.30','0.60')],
          'sample': [('S1','4.00','0.20','0.80'),('S2','5.60','0.30','0.30')]},
}
KINDS = ('blank', 'reference', 'sample')


def resolution(first, second):
    """Exact Decimal model resolution from positive Gaussian sigmas; t2 > t1."""
    t1,s1,t2,s2 = map(D, (first[1],first[2],second[1],second[2]))
    if not all(x.is_finite() for x in (t1,s1,t2,s2)) or s1<=0 or s2<=0 or t2<=t1:
        raise ValueError('Finite centers in increasing order and positive sigmas required.')
    with localcontext() as context:
        context.prec = 40
        return 2*(t2-t1)/(4*s1+4*s2)


def blank_flag(peaks):
    start,end,height = (D(RULES[k]) for k in ('blank_window_start_min','blank_window_end_min','blank_apex_flag_au'))
    return any(start<=D(p[1])<=end and D(p[3])>=height for p in peaks)


def decisions():
    rows=[]
    for case,data in CASES.items():
        first,second=data['reference']
        rs=resolution(first,second)
        passed=rs>=D(RULES['resolution_min'])
        flagged=blank_flag(data['blank'])
        rows.append({'case':case,'t_r1_min':first[1],'t_r2_min':second[1],
                     'w_b1_min':str(4*D(first[2])),'w_b2_min':str(4*D(second[2])),
                     'resolution':f'{rs:.2f}','resolution_criterion_min':RULES['resolution_min'],
                     'resolution_met':str(passed).lower(),'blank_flag':str(flagged).lower(),
                     'teaching_decision':'resolution_not_met' if not passed else 'blank_unresolved' if flagged else 'specified_checks_met',
                     'complete_sst_demonstrated':'false','origin':'synthetic'})
    return rows


def signal(peaks, time):
    return sum(float(h)*math.exp(-0.5*((time-float(mu))/float(sigma))**2) for _,mu,sigma,h in peaks)


def trace_rows():
    start,end,step=(D(RULES[k]) for k in ('trace_start_min','trace_end_min','trace_step_min'))
    count=int((end-start)/step)+1
    for case,data in CASES.items():
        for order,kind in enumerate(KINDS,1):
            for i in range(count):
                time=start+i*step
                yield [case,f'{order:02}',kind,f'{time:.2f}',f'{signal(data[kind],float(time)):.12f}','synthetic']


def write_csv(path,fields,rows):
    with path.open('w',encoding='utf-8',newline='') as stream:
        writer=csv.writer(stream,lineterminator='\n');writer.writerow(fields);writer.writerows(rows)


def write_data(folder, protect=False):
    folder=Path(folder)
    names=('sequence.csv','peak-parameters.csv','teaching-rules.csv','traces.csv','decisions.csv','simulate.py')
    if protect and any((folder/name).exists() for name in names):
        raise FileExistsError('Output exists. Choose a new output directory; existing files are preserved.')
    folder.mkdir(parents=True,exist_ok=True)
    sequence=[];parameters=[]
    for case,data in CASES.items():
        for order,kind in enumerate(KINDS,1):
            sequence.append([case,f'{order:02}',kind,len(data[kind]),'not_supplied','synthetic'])
            for label,mu,sigma,height in data[kind]:
                parameters.append([case,f'{order:02}',kind,label,mu,sigma,str(4*D(sigma)),height,'synthetic'])
    write_csv(folder/'sequence.csv',['case','injection_order','injection_kind','modeled_components','earlier_injection_history','origin'],sequence)
    write_csv(folder/'peak-parameters.csv',['case','injection_order','injection_kind','component','center_min','sigma_min','tangent_base_width_min','height_au','origin'],parameters)
    write_csv(folder/'teaching-rules.csv',['parameter','value','scope'],[[k,v,'fictional_teaching_model'] for k,v in RULES.items()])
    write_csv(folder/'traces.csv',['case','injection_order','injection_kind','time_min','signal_au','origin'],trace_rows())
    data=decisions();write_csv(folder/'decisions.csv',list(data[0]),[list(row.values()) for row in data])
    (folder/'simulate.py').write_bytes(Path(__file__).read_bytes())


if __name__=='__main__':
    parser=argparse.ArgumentParser(description=__doc__)
    parser.add_argument('--output-dir',type=Path,default=Path(__file__).resolve().parent/'reproduced',help='New output directory (default: reproduced beside this script).')
    args=parser.parse_args()
    try:write_data(args.output_dir,protect=True)
    except FileExistsError as error:parser.exit(1,str(error)+'\n')
    print('Nine synthetic injection traces and three decisions written to '+str(args.output_dir))