from pathlib import Path
import csv,datetime,importlib.util,json,math,sys
sys.path.insert(0,r'D:/Photometer_research_20261022/research/cloud_transport/private_readers')
from scipy.integrate import quad
root=Path.cwd();out=Path(r'D:/Photometer_research_20261022/revision_20261022_r4/analysis/analytic/fov_population')
out.mkdir(exist_ok=True)
spec=importlib.util.spec_from_file_location('counts_readonly',Path(__file__).with_name('count_requirements.py'))
m=importlib.util.module_from_spec(spec);spec.loader.exec_module(m)
Re=6371000.;rc=Re+17000.;rs=Re+550000.;H=rs-rc
eta=.2;dt=.003;dark=50.;width=4.;E=6.62607015e-34*299792458/337.1e-9
def geometry(theta):
 q=math.sqrt(rc*rc-rs*rs*math.sin(theta)**2)
 distance=(rs*rs-rc*rc)/(rs*math.cos(theta)+q)
 return distance,q/rc
def point_threshold(srequired,area,theta):
 distance,mu=geometry(theta)
 return srequired*4*math.pi*distance**2/(eta*area*math.cos(theta))
def footprint(full):
 half=math.radians(full/2);distance,mu=geometry(half)
 central=math.asin(rs/rc*math.sin(half))-half
 direct=2*math.pi*rc*rc*(1-math.cos(central))
 integrated=quad(lambda t:geometry(t)[0]**2/geometry(t)[1]*2*math.pi*math.sin(t),0,half,epsabs=.01)[0]
 assert abs(direct-integrated)/direct<1e-10
 return direct
rows=[];geocache={};weightcache={};flat_checks=[]
for a in [.5,2/3,1.,1.5,2.,2.5,3.]:
 for full in [10.,30.,60.,90.]:
  half=math.radians(full/2)
  # Independent closed-form flatgeometry integral.
  power=3*a-2
  closed=2*math.pi*H**(2-2*a)*(-math.log(math.cos(half)) if abs(power)<1e-12 else (1-math.cos(half)**power)/power)
  numeric=quad(lambda t:2*math.pi*H**(2-2*a)*math.sin(t)*math.cos(t)**(3*a-3),0,half,epsabs=1e-20,epsrel=1e-10)[0]
  assert abs(closed-numeric)/closed<1e-10
  flat_checks.append({'slope':a,'FOV':full,'relative_error':abs(closed-numeric)/closed})
for L in [1.7e-9,3e-9,3e-7]:
 for diameter in [1.,5.,25.]:
  area=math.pi*(diameter*.001)**2/4
  for full in [10.,30.,60.,90.]:
   half=math.radians(full/2);G=area*math.pi*math.sin(half)**2
   rate=eta*G*L*width/E
   n,s=m.count_at_power((rate+dark)*dt,m.ALPHA_OPERATIONAL,.95)
   center_threshold=point_threshold(s,area,0.)
   edge_threshold=point_threshold(s,area,half)
   surface_area=footprint(full)
   for alpha in [.5,1.,1.5,2.,2.5,3.]:
    # Dimensionless inverse source strength times surfacearea.
    weighted_area=quad(lambda t:(point_threshold(s,area,t)/1e20)**(-alpha)*geometry(t)[0]**2/geometry(t)[1]*2*math.pi*math.sin(t),0,half,epsabs=1e-4,epsrel=1e-9)[0]
    rows.append(dict(assumed_received_background_W_m2_sr_nm=L,diameter_mm=diameter,full_circular_FOV_deg=full,
     assumed_cumulative_source_number_slope_alpha=alpha,cloud_shell_altitude_km=17,orbit_altitude_km=550,
     end_to_end_efficiency=eta,projected_etendue_m2_sr=G,background_photoelectron_cps=rate,
     window_s=dt,critical_observed_count=n,required_source_count_for_95pct_detection=s,
     footprint_km2=surface_area/1e6,axis_required_isotropic_source_photons_in_window=center_threshold,
     edge_required_isotropic_source_photons_in_window=edge_threshold,
     edge_over_axis_source_threshold=edge_threshold/center_threshold,
     arbitrary_pivot_1e20ph_weighted_exposure_area_m2=weighted_area,
     hypothetical_17ns_counter_background_occupancy=rate*17e-9))
lookup={(r['assumed_received_background_W_m2_sr_nm'],r['diameter_mm'],r['full_circular_FOV_deg'],r['assumed_cumulative_source_number_slope_alpha']):r for r in rows}
for r in rows:
 L=r['assumed_received_background_W_m2_sr_nm'];d=r['diameter_mm'];f=r['full_circular_FOV_deg'];a=r['assumed_cumulative_source_number_slope_alpha']
 r['relative_detection_rate_to_same_aperture60deg']=r['arbitrary_pivot_1e20ph_weighted_exposure_area_m2']/lookup[L,d,60.,a]['arbitrary_pivot_1e20ph_weighted_exposure_area_m2']
 r['relative_detection_rate_to_25mm60deg']=r['arbitrary_pivot_1e20ph_weighted_exposure_area_m2']/lookup[L,25.,60.,a]['arbitrary_pivot_1e20ph_weighted_exposure_area_m2']
with (out/'conditional_rate_trade.csv').open('w',newline='',encoding='utf-8') as f:
 w=csv.DictWriter(f,fieldnames=list(rows[0]));w.writeheader();w.writerows(rows)
report={'recorded_utc':datetime.datetime.now(datetime.timezone.utc).isoformat(),'status':'CONDITIONAL_DESIGN_SENSITIVITY_NOT_ACTUAL_EVENT_RATE_OR_OPTIMUM',
 'questions':['B5','C8','D12','G21','G22'],
 'own_equations':{'event_law':'N_sources(>N_emitted)=C*N_emitted^(-alpha), locallyvalidbetweenallcomparedthresholds; alphaUNKNOWN,notfitted','point_counts':'S=eta*A*cos(theta)*N_emitted/(4*pi*distance^2)','surface_element':'dAcloud=distance^2/mu_cloud*dOmega_detector','diffuse_counts':'B=eta*L*DeltaLambda*A*pi*sin(FOV/2)^2*dt/Eph +dark*dt','relative_rate':'integralovercloudfootprint[N_required(theta)/pivot]^(-alpha)dAcloud; unknownC/dutyfactornotprovided','near_nadir_large_count_limit':'Ifbackgroundlimited,Srequired~sqrt(B),smallangleOmega~FOV^2, rate~diameter^alpha*Omega^(1-alpha/2); alpha2boundaryislocalapproximationonly'},
 'inputs':{'background_Ls':'Assumeduniformreceivedapparent4nmcontinuum-equivalentradiances; notmeasuredpolarclimate','eta':eta,'dark_cps':dark,'window_s':dt,'source':'Isotropicpointsourceon17kmshellwithnocloud/gasattenuation; spectralphotonnumberinacceptedband','counter':'IdealindependentlinearPoisson; highoccupancyrowsareNOTusableforecastfor17nscounter','pivot_source_photons':1e20,'pivot_meaning':'ArbitraryunitnormalizationONLY,notanassumedsourceyield; relativeanswersindependentofpivot'},
 'limits':['NoactualNBEbrightnessdistribution,geomagnetic/cloudfrequency,radio-opticalselectionoroccurrencerate','Noopticalétenduefeasibilityforchosenentrance/detector/filter; separatelyqualifiedconstraintsstillapply','FixedsourcecollectedAcos(theta),noextraunknownC(theta)maporvignetting','Nooff-axiscloudescapeanisotropy,sourceextent,windowspill,realice/aerosol/gasorcloudreflection','Countdetectiononly: usefulscienceyieldfurtherrequiresriseprecision,777upperlimit,associationandfalsepositivevalidation','Noorbitalnight/cloudfraction/attitude/safe-mode/SAA/live-timeorfollow-upcoverage; actualrateidentifiabilityremainsopen','Powerlawandspatialuniformityarehypothetical,notfitstotwelveASIMevents'],
 'hard_threshold_qualification':'Relative-ratecolumnsareBRIGHT-ENOUGHPOPULATIONWEIGHTSabovea95%-powermean-countcut, notafullintegralofP_detectoversourceluminosity. Restrictedtoabove-cutpopulation,idealexpectedtrue-detectionsliebetween.95and1timesweightup tocommonunknownnormalization; below-cutcontributionsnotincluded.',
 'validation':{'sphericalfootprint_integral_vs_analytic_cap_cases':4,'flat_closed_integrals_checked':len(flat_checks),'max_flat_relative_error':max(x['relative_error']for x in flat_checks),'rows':len(rows)},
 'examples_90_vs_60':[r for r in rows if r['diameter_mm']==25 and r['full_circular_FOV_deg']==90 and r['assumed_received_background_W_m2_sr_nm']==3e-9]}
(out/'derivation_and_qualification.json').write_text(json.dumps(report,indent=2),encoding='utf-8')
print(json.dumps({'rows':len(rows),'examples':[{'alpha':r['assumed_cumulative_source_number_slope_alpha'],'relative_90_to_60':r['relative_detection_rate_to_same_aperture60deg'],'edge_axis':r['edge_over_axis_source_threshold']}for r in report['examples_90_vs_60']],'validation':report['validation']},indent=2))
