from atmPy.aerosols.instruments.POPS import mie as pm

single_pnt_cali_d = 500
single_pnt_cali_d *= 1e-3
d,mie = pm.makeMie_diameter(diameterRangeInMikroMeter=[single_pnt_cali_d,single_pnt_cali_d + 1e-3],noOfdiameters=2)

d,mie

%matplotlib inline

reload(calibration)

cal,a  = calibration.generate_calibration(
    single_pnt_cali_d=508,
    single_pnt_cali_ior=1.6,
    single_pnt_cali_int=1000,
    ior=complex(1.5,0.01),
    dr=[100, 5000],
    no_pts=500,
    no_cal_pts=40,
    plot=True,
    raise_error=True,
    test=False,
)

single_pnt_cali_d = 500

mie_cal = cal.data

dt = mie_cal.index[abs(mie_cal.index - single_pnt_cali_d).argmin()]
dt

at =mie_cal.loc[dt,'amp']


diter = cal.calibrationFunction(at)
step = 0.01


from scipy import optimize

amp = optimize.fsolve(lambda x: cal.calibrationFunction(x) - single_pnt_cali_d, at)

cal.calibrationFunction(amp)

mie_cal

calold.calibrationFunction(500)

dft  = calibration.generate_calibration(
    single_pnt_cali_d=508,
    single_pnt_cali_ior=1.6,
    single_pnt_cali_int=1000,
    ior=2.95,
    dr=[100, 5000],
    no_pts=500,
    no_cal_pts=40,
    plot=True,
    raise_error=True,
    test=False,
)

dft

dr=[100, 5000]
drum = np.array(dr)/1e3
d, amp = pm.makeMie_diameter(diameterRangeInMikroMeter=drum,
                                 noOfdiameters=500,
                                 IOR=2.95)

plt.plot(d,amp)
plt.loglog()

cal  = calibration.generate_calibration(
    single_pnt_cali_d=508,
    single_pnt_cali_ior=1.6,
    single_pnt_cali_int=1000,
    ior=1.5,
    dr=[100, 5000],
    no_pts=500,
    no_cal_pts=40,
    plot=False,
    raise_error=True,
    test=False,
)

cal.plot()

np.round()

cal.plot()

ns = np.linspace(1.5,2.95,3)
ns = [2.]

f,a = plt.subplots()
a.loglog()
for n in ns:
    d,mie = pm.makeMie_diameter(IOR = n)
    g, = a.plot(d,mie, label = n)
#     g.set_linestyle('')
#     g.set_marker('.')
a.legend()

%%time
n = 2.95
d,amp = pm.makeMie_diameter(diameterRangeInMikroMeter=[0.1, 5.0],  
                            noOfdiameters=1e4,
                            IOR = n)

df = pd.DataFrame({'d':d, 'amp':amp})

# %%time
nofbins = 50
valid = False
while not valid:
    nofbins -= 1
    binedgs = np.logspace(np.log10(mie.min()), np.log10(mie.max()), nofbins)
    binedgs[0] -= (binedgs[1] - binedgs[0]) * 0.01 # this is to ensure the first point is not onthe edge ... for cut function used later


    dfavg = df.groupby(pd.cut(df.amp, binedgs)).median()
    dfstd = df.groupby(pd.cut(df.amp, binedgs)).mad()
    
    dfavg.index = dfavg.d
    dfstd.index = dfavg.d
    dfavg['sigma_d'] = dfstd.d
    
    dfavg.index.name = None
    # dfbined

    dfavg.sort_values('amp', axis=0, inplace=True)
    if not ((dfavg.d.values[1:] - dfavg.d.values[:-1]) < 0).sum():
        valid = True
if 0:
    a = df.mie.hist(bins = dfavg)
    a.semilogx()
nofbins

f,a = plt.subplots()
g, = a.plot(d,amp, label = n)
g.set_linestyle('')
g.set_marker('.')
g.set_markersize(1)

#print binedges
if 1:
    for be in binedgs:
        out = a.axhline(be)
        out.set_linestyle(':')
        out.set_linewidth(1)
        out.set_color('0.5')
        
#print averaged data
if 1:
#     a.plot(dfavg.d, dfavg.mie)
#     a.errorbar(dfavg.d, dfavg.amp, xerr=dfstd.d, yerr=dfstd.amp, zorder = 10, ls ='')
    a.errorbar(dfavg.d, dfavg.amp, 
               xerr=dfavg.sigma_d,
#                xerr = dfstd.amp,
#                yerr=dfstd.amp, 
               zorder = 1, ls ='')
    
#test
if 1:
    plt.plot(dnew, ampnew)

#     plt.plot(d,amp)


a.legend()
a.loglog()


dfstd.amp

dfavg.sigma_d

reload(calibration)

# %%debug -b /mnt/telg/programs/anaconda3/lib/python3.7/site-packages/scipy/interpolate/fitpack2.py:177
cali = calibration.Calibration(dfavg)

cali.plot_calibration()

bla = type("bla", (), {})
self = bla()
self.data = dfavg

from scipy import interpolate

fkt = interpolate.interp1d(self.data.amp, self.data.d, kind = 'quadratic')

ampnew = np.logspace(np.log10(self.data.amp.min()), np.log10(self.data.amp.max()), 100)
dnew = fkt(ampnew)

plt.plot(self.data.d, self.data.amp)
plt.plot(dnew, ampnew)
plt.loglog()

sf = UnivariateSpline(np.log10(self.data.d.values), np.log10(self.data.amp.values), 
                      s=0,
                      k = 2
                     )

d = np.linspace(np.log10(self.data.d.values.min()), np.log10(self.data.d.values.max()), 500)
amp = sf(d)

plt.plot(d,amp)
# plt.plot(np.log10(dfavg.d), np.log10(dfavg.amp))

# us = UnivariateSpline(np.log10(self.data.amp), np.log10(self.data.d), s=0)
us = UnivariateSpline(amp, d, s=0)
cal_function = lambda amp: 10**us(np.log10(amp))

from scipy.interpolate import UnivariateSpline

UnivariateSpline()

df.shape bide

df.mie.mod(binedgs)

"""
Explanation: (baustelle) Generate a theoretical calibration
End of explanation
"""


cali, a = calibration.generate_calibration(single_pnt_cali_d = 510,
                                             single_pnt_cali_ior = 1.6,
                                             single_pnt_cali_int = 10**3.508,
                                             ior = 1.7,
                                             dr=np.array([110, 3400]),
                                             no_pts = 600,
                                             no_cal_pts = 29,
                                             plot = True,
                                             raise_error = True,
                                             test = False
                                            )

cali.plot_calibration()

"""
Explanation: One can generate a calibration file based on the theoretical responds of the POPS instrument and a single point calibration.
End of explanation
"""


cali, a = calibration.generate_calibration(single_pnt_cali_d = 510,
                                             single_pnt_cali_ior = 1.6,
                                             single_pnt_cali_int = 10**3.508,
                                             ior = 1.5,
                                             dr=np.array([110, 3400]),
                                             no_pts = 600,
                                             no_cal_pts = 30,
                                             plot = True,
                                             raise_error = True,
                                             test = False
                                            )

cali.plot_calibration()

"""
Explanation: Generate a theoretical calibration
One can generate a calibration file based on the theoretical responds of the POPS instrument and a single point calibration.
End of explanation
"""


data = """130 1.73
140 1.85
150 1.99
173 2.28
200 2.53
233 2.77
270 2.96
315 3.04
365 3.165
420 3.213
490 3.355
570 3.56
660 3.72
770 3.795
890 3.86
1040 4.0
1200 4.05
1400 4.1
1600 4.18
1880 4.32
2180 4.44
2500 4.5"""

cal = calibration.read_str(data)

"""
Explanation: Generate an experimental calibration
This section shows how to read and save calibration data. In the process the data will be fit and a calibration function is created. Data and function can be plotted using the plot function mentioned below.
read from string
create a string where each line consist of a diameter value (in nm) and an amplitude (in digitizer bin) seperated by space. <br> Note, in the example below the secend column is the $log_{10}(amplitude)$. The following cell does the convertion!
End of explanation
"""


fname = './data/POPS_calibration.csv'
cal.save_csv(fname)

"""
Explanation: Save the calibration from above
End of explanation
"""


out = calibration.read_csv(fname)

"""
Explanation: read the file we just saved
End of explanation
"""


f,a,df,gf = cal.plot_calibration()

"""
Explanation: plot the created calibration function and data
End of explanation
"""


out = calibration.get_interface_bins('/Users/htelg/data/POPS_calibrations/150529_MS1_DOS.csv', 16)

"""
Explanation: Translate bins to diameters
This information is necessarey if data is provided prebinned e.g. when recorded through the serial interface.
End of explanation
"""


