"""Independent public discrete native-node population, with numerical cell integrals."""
import numpy as np
from scipy.integrate import quad
def produce(native,seed=32032,slope=2.3,binary=.35,concentration=.4,aperture=.3,distance=850):
    rows=np.array(native['rows']);nodes=rows[:,0];edges=np.r_[nodes[0],(nodes[:-1]+nodes[1:])/2,nodes[-1]]
    norm=quad(lambda mass:mass**(-slope),nodes[0],nodes[-1],epsabs=1e-13)[0]
    weights=np.array([quad(lambda mass:mass**(-slope),lo,hi,epsabs=1e-13)[0]/norm for lo,hi in zip(edges[:-1],edges[1:])]);cdf=np.cumsum(weights);cdf[-1]=1
    state=int(seed)
    def random():
        nonlocal state
        state=(1664525*state+1013904223)%2**32
        return (state+.5)/2**32
    node_ids=np.array([np.searchsorted(cdf,random(),side='right') for _ in range(1024)])
    masses=nodes[node_ids];ids=list(range(1024))
    for index in range(1023,0,-1):
        target=int(random()*(index+1));ids[index],ids[target]=ids[target],ids[index]
    pairs=int(np.floor(binary*1024/(1+binary)));systems=[ids[2*i:2*i+2] for i in range(pairs)]+[[item] for item in ids[2*pairs:]]
    position_latents=[]
    for _ in range(1024):
        u,v,w=random(),random(),random();radius=(u**(-2/3)-1)**(-.5);theta=np.arccos(2*v-1);phi=2*np.pi*w
        position_latents.append(radius*np.array([np.sin(theta)*np.cos(phi),np.sin(theta)*np.sin(phi),np.cos(theta)]))
    records=[];all_light=0;aperture_light=0
    for index,components in enumerate(systems):
        mass=float(sum(masses[components]));position=position_latents[components[0]]*2*(mass/.5)**(-concentration);ray=position+np.array([0,0,distance]);actual_distance=np.linalg.norm(ray)
        angle=np.degrees(np.arctan2(np.hypot(ray[0],ray[1]),ray[2]));magnitudes=rows[node_ids[components],5:8]
        flux=np.exp(-.4*np.log(10)*magnitudes).sum(axis=0);apparent=-2.5/np.log(10)*np.log(flux)+5*np.log10(actual_distance/10);in_aperture=bool(angle<aperture)
        light=flux[0]*(10/actual_distance)**2;all_light+=light
        if in_aperture:aperture_light+=light
        records.append({'id':index,'componentIds':components,'initialMasses':masses[components].tolist(),'initialSystemMass':mass,'positionPc':position.tolist(),'distancePc':float(actual_distance),'angularRadiusDegrees':float(angle),'componentMagnitudes':magnitudes.tolist(),'apparentMagnitudes':apparent.tolist(),'inAperture':in_aperture,'detected':bool(in_aperture and apparent[0]<18.5)})
    return {'schema':'public-node-reference-v1','systems':records,'componentMasses':masses.tolist(),'nodeWeights':weights.tolist(),'supportedNativeAgeYears':native['ageYears'],'passbands':native['passbands'],'apertureLightFraction':float(aperture_light/all_light),'discreteBirthMass':float(masses.sum()),'continuousIMFMeanMass':quad(lambda mass:mass*mass**(-slope),nodes[0],nodes[-1])[0]/norm,'discreteExpectedMeanMass':float(np.dot(weights,nodes)),'computedDynamicalEvolution':False}