"""Independent positive time quadrature over saved fine masses; no cohort-mean interpolation in log space."""
import math
import numpy as np
from scipy.special import logsumexp
from evolution import histories

def nodes_from_saved_masses(record,fixed,order):
 time=np.array(record['timeGyr']);states=np.array(record['states']);nodes,weights=np.polynomial.legendre.leggauss(order);dt=np.diff(time);fraction=(nodes+1)/2
 clock=(time[:-1,None]+dt[:,None]*fraction).reshape(-1);mass=(states[:-1,None,:]+(states[1:]-states[:-1])[:,None,:]*fraction[None,:,None]).reshape(-1,14)
 mg=mass[:,4]/mass[:,0];fe=mass[:,5]/mass[:,0]
 if np.any(mg<=0)or np.any(fe<=0):raise ValueError('Positive earlyMg/Fe support lost')
 means=np.column_stack((np.log10(fe/fixed['solarMassFractions']['fe']),np.log10((mg/fe)/(fixed['solarMassFractions']['mg']/fixed['solarMassFractions']['fe']))))
 weighted=(dt[:,None]*weights[None,:]/2).reshape(-1)*np.array([histories(t,record['parameters'])[0]for t in clock]);normalization=float(weighted.sum());weighted/=normalization
 return means,weighted,normalization

def quadrature_logdensity(record,fixed,rows,rho,order,fe_shift=0,ratio_shift=0,prepared=None):
 means,weights,total=prepared if prepared is not None else nodes_from_saved_masses(record,fixed,order);means=means+np.array([fe_shift,ratio_shift]);answer=[]
 for row in rows:
  error=np.array([row['fehError'],row['mgfeError']]);r=(np.array([row['feh'],row['mgfe']])-means)/error
  value=-.5*(r[:,0]**2-2*rho*r[:,0]*r[:,1]+r[:,1]**2)/(1-rho*rho)-math.log(2*math.pi*error.prod()*math.sqrt(1-rho*rho));answer.append(float(logsumexp(value+np.log(weights))))
 return np.array(answer),total