export const MODEL_VERSION='1.0.0';
export const DEFAULT_PARAMETERS=Object.freeze({period:1.63,beta:.12,forcing:.02,ratio:1,initialAngle:2*Math.PI/180});
export const PARAMETER_BOUNDS=Object.freeze({period:[1.5,1.8],beta:[.05,.4],forcing:[.005,.02],ratio:[.8,1.2],initialAngle:[0,2*Math.PI/180]});
export function validateParameters(p){for(const[k,[lo,hi]]of Object.entries(PARAMETER_BOUNDS))if(!Number.isFinite(p[k])||p[k]<lo-1e-14||p[k]>hi+1e-14)throw new RangeError(`${k} outside [${lo},${hi}]`);return p;}
export function response(p){validateParameters(p);const omega0=2*Math.PI/p.period,omega=omega0*p.ratio,real=omega0*omega0-omega*omega,imag=2*p.beta*omega;return{omega0,omega,omegaD:Math.sqrt(omega0*omega0-p.beta*p.beta),amplitude:p.forcing/Math.hypot(real,imag),phase:Math.atan2(-imag,real),peakRatio:Math.sqrt(1-2*(p.beta/omega0)**2)};}
export function motion(p,t,mode='forced'){const q=response(p);if(!Number.isFinite(t)||t<0||!['free','forced'].includes(mode))throw new RangeError('invalid time/mode');const A=mode==='forced'?q.amplitude:0,phase=q.phase,C=p.initialAngle-A*Math.cos(phase),D=(p.beta*C+A*q.omega*Math.sin(phase))/q.omegaD,cos=Math.cos(q.omegaD*t),sin=Math.sin(q.omegaD*t),decay=Math.exp(-p.beta*t),theta=A*Math.cos(q.omega*t+phase)+decay*(C*cos+D*sin),velocity=-A*q.omega*Math.sin(q.omega*t+phase)+decay*((-p.beta*C+q.omegaD*D)*cos+(-p.beta*D-q.omegaD*C)*sin),drive=mode==='forced'?p.forcing*Math.cos(q.omega*t):0;return{theta,velocity,drive,acceleration:drive-2*p.beta*velocity-q.omega0*q.omega0*theta,energy:.5*(velocity*velocity+q.omega0*q.omega0*theta*theta)};}
export function decayPeaks(p){const q=response(p);if(p.initialAngle<=0)throw new RangeError('自由衰减需正的初始角');return Array.from({length:10},(_,n)=>{const time=2*Math.PI*n/q.omegaD;return{time,amplitude:p.initialAngle*Math.exp(-p.beta*time)};});}
export function measurement(p,kind,raw,source='ideal-model'){validateParameters(p);if(!['decay','response'].includes(kind)||!['ideal-model','manual-reading'].includes(source))throw new RangeError('unknown kind/source');if(!Number.isFinite(raw.amplitude)||raw.amplitude<=0||raw.amplitude>10*Math.PI/180)throw new RangeError('幅度须为0–10°之间的正数');if(kind==='decay'&&(!Number.isFinite(raw.time)||raw.time<0))throw new RangeError('峰时刻须为非负数');if(kind==='response'&&(!Number.isFinite(raw.phase)||raw.phase<-Math.PI||raw.phase>0))throw new RangeError('相位须在−180°至0°内');return{experiment:'course-02/pohl-resonance',modelVersion:MODEL_VERSION,kind,source,parameters:{...p},raw:{...raw},derived:kind==='decay'?{logAmplitude:Math.log(raw.amplitude)}:{frequency:response(p).omega/(2*Math.PI),ratio:p.ratio}};}
export function sweep(p){return Array.from({length:21},(_,i)=>{const q={...p,ratio:.8+i*.02},r=response(q);return measurement(q,'response',{amplitude:r.amplitude,phase:r.phase});});}
const same=(a,b,exclude=[])=>Object.keys(PARAMETER_BOUNDS).filter(k=>!exclude.includes(k)).every(k=>a[k]===b[k]);
export function fitDecay(rows,p,source){const q=rows.filter(r=>r.kind==='decay'&&r.source===source&&same(r.parameters,p)),times=[...new Set(q.map(r=>r.raw.time))];if(times.length<3)throw new RangeError('至少记录三个不同时刻的同条件峰');const x=q.map(r=>r.raw.time),y=q.map(r=>Math.log(r.raw.amplitude)),xm=x.reduce((a,b)=>a+b,0)/x.length,ym=y.reduce((a,b)=>a+b,0)/y.length,sxx=x.reduce((a,b)=>a+(b-xm)**2,0),slope=x.reduce((a,b,i)=>a+(b-xm)*(y[i]-ym),0)/sxx,intercept=ym-slope*xm;return{beta:-slope,intercept,slope,x,y,count:q.length};}
export function responseRecords(rows,p,source){return rows.filter(r=>r.kind==='response'&&r.source===source&&same(r.parameters,p,['ratio'])).sort((a,b)=>a.parameters.ratio-b.parameters.ratio);}