Add Support for a calculation of LCL

[xebadir]

For the SB parcel one can just use the Stull (2000) approximation: SBLCL= 125.*(((T[0,:,:]+T[1,:,:])/2)-Td[0,:,:])

For the ML case - you need to mix the parcel (e.g. take the mean of the layer characteristics) MLLCL_50mb = np.zeros((ny,nx)) MLLCL_100mb = np.zeros((ny,nx)) for i in range(0,nx): for j in range(0,ny): P_in = P[:,j,i] T_in = T[:,j,i] Td_in = Td[:,j,i] MU = int(MUlvl[j,i]) MLLCL_50mb[j,i],Q_ML50[j,i],avgt50[j,i],avgtd50[j,i],avgq50[j,i],avgp50[j,i]= MLLCL.getmllcl(p_in=P_in,t_in=T_in,td_in=Td_in,ml_depth=1000.0,nk=nz) MLLCL_100mb[j,i],Q_ML100[j,i],avgt100[j,i],avgtd100[j,i],avgq100[j,i],avgp100[j,i]= MLLCL.getmllcl(p_in=P_in,t_in=T_in,td_in=Td_in,ml_depth=1000.0,nk=nz) Q_MU[j,i]=Q[MU,j,i] H_MU[j,i]=H[MU,j,i] P_MU[j,i]=P[MU,j,i] Though this was for consistency and pure convenience using a stripped back version of George's code. In any case, you calculate the LCL using the same calculation as above, just for the mixed properties.

There is also an analytical python/fortran solution by David Romps:

surfRH = calc.relative_humidity_from_dewpoint(T[0,:,:]*units.units('degC'),
                 Td[0,:,:]*units.units('degC')).magnitude
        ML50_RH = calc.relative_humidity_from_dewpoint(avgt50[:,:]*units.units('degK'),
                 avgtd50[:,:]*units.units('degK')).magnitude
        ML100_RH = calc.relative_humidity_from_dewpoint(avgt100[:,:]*units.units('degK'),
                 avgtd100[:,:]*units.units('degK')).magnitude
        #t0 = time.clock() 
        #for i in range(0,nx):
        #   for j in range(0,ny):
        #       exact_SBLCL[j,i] = lcl(P[0,j,i]*100,T[0,j,i]+273.15,rh=surfRH[j,i]) 
        #print time.clock() - t0, "seconds full grid SBLCL Romps Py Calc time"
        Pin = P[0,:,:]*100
        Tin = T[0,:,:]+273.15
        t0 = time.clock()
        for i in range(0,nx):
           for j in range(0,ny):
                exact_SBLCL[j,i] = exactLCL.lcl_mod.lcl(Pin[j,i],Tin[j,i],surfRH[j,i])
                exact_MLLCL50[j,i] = exactLCL.lcl_mod.lcl(avgp50[j,i],avgt50[j,i],ML50_RH[j,i])
                exact_MLLCL100[j,i] = exactLCL.lcl_mod.lcl(avgp100[j,i],avgt100[j,i],ML100_RH[j,i])

I'll try to implement this - but the Stull approximate will do in a pinch.

[chiaral]

For the ML case - you need to mix the parcel (e.g. take the mean of the layer characteristics) MLLCL_50mb = np.zeros((ny,nx)) MLLCL_100mb = np.zeros((ny,nx)) for i in range(0,nx): for j in range(0,ny): P_in = P[:,j,i] T_in = T[:,j,i] Td_in = Td[:,j,i] MU = int(MUlvl[j,i]) MLLCL_50mb[j,i],Q_ML50[j,i],avgt50[j,i],avgtd50[j,i],avgq50[j,i],avgp50[j,i]= MLLCL.getmllcl(p_in=P_in,t_in=T_in,td_in=Td_in,ml_depth=1000.0,nk=nz) MLLCL_100mb[j,i],Q_ML100[j,i],avgt100[j,i],avgtd100[j,i],avgq100[j,i],avgp100[j,i]= MLLCL.getmllcl(p_in=P_in,t_in=T_in,td_in=Td_in,ml_depth=1000.0,nk=nz) Q_MU[j,i]=Q[MU,j,i] H_MU[j,i]=H[MU,j,i] P_MU[j,i]=P[MU,j,i] Though this was for consistency and pure convenience using a stripped back version of George's code. In any case, you calculate the LCL using the same calculation as above, just for the mixed properties.

I have this code. It should be easy to wrap it with my loops for model lev, pressure lev, and pressure lev 1d and then keep wrapping it in the same structure.

[chiaral]

The problem will be to define tests. do we have ways to test this outside of our codes? metpy?

[chiaral]

Q_MU[j,i]=Q[MU,j,i] H_MU[j,i]=H[MU,j,i] P_MU[j,i]=P[MU,j,i]

Does it need to be H or AGLH? right now i simply calculate AGLH, by fixing the surface height to 2m

[xebadir]

Should be AGL heights.

[xebadir]

Should be AGL heights.