Review index · Download raw source

reference/statistics/browser_math.js

1(function(root){'use strict';2const H=6.62607015e-34,C=299792458,ER=H*C/(337.1e-9),RPH=1e10/(4*Math.PI);3function lgamma(z){if(!Number.isFinite(z)||z<=0)throw Error('Log-gamma requires a positive finite argument');const p=[676.5203681218851,-1259.1392167224028,771.32342877765313,-176.61502916214059,12.507343278686905,-.13857109526572012,9.9843695780195716e-6,1.5056327351493116e-7];if(z<.5)return Math.log(Math.PI)-Math.log(Math.sin(Math.PI*z))-lgamma(1-z);z-=1;let x=.99999999999980993;for(let i=0;i<p.length;i++)x+=p[i]/(z+i+1);let t=z+7.5;return .9189385332046727+(z+.5)*Math.log(t)-t+Math.log(x)}4function mass(k,mu){if(!Number.isFinite(k)||!Number.isInteger(k)||k<0||!Number.isFinite(mu)||mu<0)throw Error('Invalid Poisson mass arguments');if(mu===0)return k===0?1:0;if(k===0)return Math.exp(-mu);return Math.exp(-mu+k*Math.log(mu)-lgamma(k+1))}5function sf(n,mu){if(!Number.isFinite(n)||!Number.isFinite(mu)||mu<0)throw Error('Invalid Poisson tail arguments');n=Math.ceil(n);if(n<=0)return 1;if(mu===0)return 0;if(n===1)return -Math.expm1(-mu);if(mu<0||!Number.isFinite(mu))throw Error('Poisson mean must be finite and nonnegative');let k,term,s;if(n-1<mu){k=n-1;term=mass(k,mu);s=term;while(k>0){term*=k/mu;k--;s+=term;if(term<1e-15*s)break}return Math.max(0,Math.min(1,1-s))}k=n;term=mass(k,mu);s=term;for(let j=0;j<200000;j++){k++;term*=mu/k;s+=term;if(term<1e-15*s||term===0)break;if(j===199999)throw Error('Tail did not converge')}return Math.max(0,Math.min(1,s))}6function cdf(n,mu){if(!Number.isFinite(n)||!Number.isFinite(mu)||mu<0)throw Error('Invalid Poisson CDF arguments');if(n<0)return 0;if(mu===0)return 1;n=Math.floor(n);if(n===0)return Math.exp(-mu);if(n<mu){let k=n,t=mass(k,mu),s=t;while(k>0){t*=k/mu;k--;s+=t;if(t<1e-15*s)break}return Math.max(0,Math.min(1,s))}return 1-sf(n+1,mu)}7function critical(B,alpha){if(!Number.isFinite(B)||B<0)throw Error('Invalid background mean');if(!Number.isFinite(alpha)||!(alpha>0&&alpha<.5))throw Error('Invalid false-alarm allocation');if(B>1e6)throw Error('Exact live threshold is limited to background means ≤1,000,000; use the offline algorithm for larger cases');let lo=0,hi=Math.max(1,Math.ceil(B+12*Math.sqrt(B+1)+60));while(sf(hi,B)>alpha)hi*=2;while(lo+1<hi){const m=Math.floor((lo+hi)/2);if(sf(m,B)<=alpha)hi=m;else lo=m}return hi}8function threshold(B,alpha,power=.95){if(!Number.isFinite(power)||power<=0||power>=1)throw Error('Detection power must lie between zero and one');const n=critical(B,alpha);let lo=B,hi=Math.max(B+10,n+10*Math.sqrt(n+1));while(sf(n,hi)<power)hi*=2;for(let i=0;i<65;i++){const m=(lo+hi)/2;if(sf(n,m)<power)lo=m;else hi=m}return{critical:n,signal:Math.max(0,hi-B),power:sf(n,hi),tail:sf(n,B)}}9function upper(n,B=0,alpha=.05){if(!Number.isInteger(n)||!Number.isFinite(B)||B<0||!Number.isFinite(alpha)||alpha<=0||alpha>=1)throw Error('Upper limit requires integer counts, finite background and0<alpha<1');if(n>1e6||B>1e6)throw Error('Use offline upper-limit algorithm above one million counts');n=Math.floor(n);if(n<0||B<0)throw Error('Counts and background must be nonnegative');let lo=0,hi=Math.max(10,n+10*Math.sqrt(n+1));while(cdf(n,hi)>alpha)hi*=2;for(let i=0;i<70;i++){let m=(lo+hi)/2;if(cdf(n,m)>alpha)lo=m;else hi=m}return Math.max(0,hi-B)}10function budget({d=5,fov=60,eff=.2,band=4,dt=.003,L=1.77e-9,aurora=0,fill=1,dark=50,fluence=1e8}){if(![d,fov,eff,band,dt,L,aurora,fill,dark,fluence].every(Number.isFinite)||!(d>0&&fov>0&&fov<180&&eff>0&&eff<=1&&band>0&&dt>0&&L>=0&&aurora>=0&&fill>=0&&fill<=1&&dark>=0&&fluence>=0))throw Error('Inputs outside physical bounds');const A=Math.PI*(d*.001)**2/4,G=A*Math.PI*Math.sin(fov*Math.PI/360)**2,sky=eff*G*(L*band/ER+aurora*fill*RPH),B=(sky+dark)*dt,S=eff*A*fluence;if(![A,G,sky,B,S].every(Number.isFinite))throw Error('Input magnitudes overflow the numerical budget domain');return{A,G,sky,B,S,energy:ER,occupancy_background_17ns:(sky+dark)*17e-9,shot_sigma:Math.sqrt(B)}}11function filterCenter(angle,neff=1.8,wave=337.1){if(![angle,neff,wave].every(Number.isFinite)||angle<0||angle>=90||neff<=1||wave<=0)throw Error('Filter model requires0<=AOI<90, effective index>1 and positive wavelength');return wave*Math.sqrt(1-Math.sin(angle*Math.PI/180)**2/neff**2)}12function capture(sw,width,angle,neff=1.8){if(!Number.isFinite(sw)||!Number.isFinite(width)||sw<=0||width<=0)throw Error('Source and filter widths must be positive and finite');let ss=sw/2.3548200450309493,sf=width/2.3548200450309493,delta=filterCenter(angle,neff)-337.1;return sf/Math.sqrt(ss*ss+sf*sf)*Math.exp(-delta*delta/(2*(ss*ss+sf*sf)))}13const M={H,C,ER,RPH,lgamma,sf,cdf,critical,threshold,upper,budget,filterCenter,capture};root.PhotometerMath=M;if(typeof module!=='undefined'&&module.exports)module.exports=M;14})(typeof globalThis!=='undefined'?globalThis:this);15