In [1]:
from netCDF4 import Dataset
import numpy as np
import sys
import time

In [2]:
flt = ['20150617']

flts = ['20150617','20150620','20150701','20150702','20150706','20150709']

savePath = '/Users/danstechman/GoogleDrive/PECAN-Data/'

# Recreate all parameter files
runAll = True

In [3]:
if runAll:
    runFlights = flts
else:
    runFlights = flt

In [4]:
for flight in runFlights:
    if flight == '20150617':
        ## Flight-level data parameters
        flFile = '20150617I1_AXC.nc'
        dewPtSens = 'TDM.1'
        altSens = 'AltGPS.2'
        wsSens = 'WS.dX'


        # Start and end times of spirals in seconds
        startT = np.array([14203, 16750, 18576, 20353, 21846, 23588, 24995],dtype=float)
        endT = np.array([14917, 17948, 19577, 21440, 22905, 24570, 25895],dtype=float)
        
        # Start and end times of PDDs in seconds
        PDDstartT = np.array([9272, 9671, 10827, 12996, 14925, 16023, 19682, 21499, 24672, 26046],dtype=float)
        PDDendT = np.array([9647, 10433, 11843, 14107, 15484, 16602, 20322, 21846, 24995, 26376],dtype=float)

        # Time/[temp at same time] when last evidence of ice was observed
        mlBotTemp = np.array([np.nan, 2.3943, np.nan, 3.3489, 4.3656, 2.5243, 5.2405],dtype=float)
        mlBotTime = np.array([np.nan, 17018, np.nan, 20530, 22632, 23810, 25750],dtype=float)
        # S1 - No particles warmer than -5
        # S2 - No particles between +2 and +8, so ML bottom actually warmer than +2
        # S3 - No particles warmer than -4
        
        # Time/[temp at same time] when first evidence of liquid water was observed
        mlTopTemp = np.array([np.nan, 0.3578, np.nan, 0.9657, 0.4826, 0.1111, 0.3433],dtype=float)
        mlTopTime = np.array([np.nan, 17074, np.nan, 20580, 22512, 23857, 25635],dtype=float)
        # S1 - No particles warmer than -5
        # S3 - No particles warmer than -4

        # MCS type
        # Defined based on general MCS type at time of each spiral
        # T = Trailing Stratiform, L = Leading Stratiform,
        # P = Parallel Stratiform, C = Cluster MCS, F = Post-Frontal
        mcsType = np.array(['T','T','T','T','T','T','T'],dtype=str)

        # Location of spiral relative to MCS structure. 
        # T = Transition Zone, S = Stratiform Region, A = Anvil Region, X = N/A
        sprlZone = np.array(['A', 'A', 'A', 'S', 'S', 'S', 'S'],dtype=str)

        # Thresholds of particle interarrival time (in seconds) for each spiral 
        # below which particles are to be considered shattered
        CIP_intArvThrsh = np.array([1.9e-5, 1.2e-5, 7.8e-5, 3.7e-6, 7e-6, 8.7e-6, 3e-6],dtype=float)
        
        PIP_rjctStartT = np.array([np.nan, np.nan, np.nan, np.nan, np.nan, np.nan, np.nan],dtype=float)
        PIP_rjctEndT = np.array([np.nan, np.nan, np.nan, np.nan, np.nan, np.nan, np.nan],dtype=float)

    elif flight == '20150620':
        flFile = '20150620I1_AXC.nc'
        dewPtSens = 'TDM.1'
        altSens = 'AltGPS.2'
        wsSens = 'WS.dX'

        startT = np.array([17757, 18230, 19284, 20845, 21877, 23890, 27947],dtype=float)
        endT = np.array([18229, 18783, 20227, 21333, 22968, 24451, 28702],dtype=float)
        
        PDDstartT = np.array([15109, 17115, 18841, 20289, 22907, 24433, 25740, 26898],dtype=float)
        PDDendT = np.array([17002, 17715, 19295, 20838, 23876, 25709, 26774, 27973],dtype=float)

        mlBotTemp = np.array([np.nan, 6.8400, np.nan, 5.9069, 5.7851, 4.6184, 4.7246],dtype=float)
        mlBotTime = np.array([np.nan, 18647, np.nan, 21219, 22080, 24352, 28124],dtype=float)
        # S1 - No particles warmer than -9
        # S2 - Very little melting all the way through +7 (nothing recorded warmer)
        # S3 - No particles warmer than -6
        
        mlTopTemp = np.array([np.nan, 0.0282, np.nan, 0.1974, 0.1935, 2.2460, 0.0326],dtype=float)
        mlTopTime = np.array([np.nan, 18582, np.nan, 21122, 22210, 24276, 28268],dtype=float)
        # S1 - No particles warmer than -9
        # S3 - No particles warmer than -6

        mcsType = np.array(['T','T','T','T','T','T','T'],dtype=str)

        sprlZone = np.array(['T', 'T', 'A', 'S', 'S', 'S', 'S'],dtype=str)

        CIP_intArvThrsh = np.array([1.5e-5, 7e-6, 3e-5, 6e-6, 6e-6, 1.5e-5, 3e-6],dtype=float)
        
        PIP_rjctStartT = np.array([np.nan, np.nan, np.nan, np.nan, np.nan, np.nan, np.nan],dtype=float)
        PIP_rjctEndT = np.array([np.nan, np.nan, np.nan, np.nan, np.nan, np.nan, np.nan],dtype=float)
        
    elif flight == '20150701':
        flFile = '20150701I1_AXC.nc'
        dewPtSens = 'TDM.2X'
        altSens = 'AltGPS.2'
        wsSens = 'WS.d'

        startT = np.array([21882],dtype=float)
        endT = np.array([22838],dtype=float)
        
        PDDstartT = np.array([17421, 21247, 22838, 23775, 25230, 25899, 26398, 27598, 28944, 30158, 31008],dtype=float)
        PDDendT = np.array([20037, 21901, 23699, 25152, 25819, 26318, 27519, 28871, 30052, 30978, 31920],dtype=float)

        mlBotTemp = np.array([3.4335],dtype=float)
        mlBotTime = np.array([22482],dtype=float)
        
        mlTopTemp = np.array([0.7108],dtype=float)
        mlTopTime = np.array([22375],dtype=float)

        mcsType = np.array(['C'],dtype=str)

        sprlZone = np.array(['S'],dtype=str)

        CIP_intArvThrsh = np.array([3e-5],dtype=float)
        
        PIP_rjctStartT = np.array([np.nan],dtype=float)
        PIP_rjctEndT = np.array([np.nan],dtype=float)
        
    elif flight == '20150702':
        flFile = '20150702I1_AXC.nc'
        dewPtSens = 'TDM.2X'
        altSens = 'AltGPS.3'
        wsSens = 'WS.d'

        startT = np.array([13336, 15350, 16892],dtype=float)
        endT = np.array([14551, 16501, 18196],dtype=float)
        
        PDDstartT = np.array([10203, 12967, 14692, 16500, 18177, 20573, 22860],dtype=float)
        PDDendT = np.array([11734, 13335, 15350, 16932, 20497, 22783, 25177],dtype=float)

        mlBotTemp = np.array([np.nan, 1.5034, np.nan],dtype=float)
        mlBotTime = np.array([np.nan, 15826, np.nan],dtype=float)
        # S1 - No particles warmer than -13
        # S3 - No particles between ~+1 and +9.59 - ML bottom somewhere in here
        
        mlTopTemp = np.array([np.nan, 0.7152, np.nan],dtype=float)
        mlTopTime = np.array([np.nan, 15852, np.nan],dtype=float)
        # S1 - No particles warmer than -13
        # S3 - Melting not observed at T<0.79; no particles warmer than this before +9.59

        mcsType = np.array(['L','L','L'],dtype=str)

        sprlZone = np.array(['A', 'A', 'A'],dtype=str)

        CIP_intArvThrsh = np.array([1.7e-4, 3.9e-5, 3.9e-5],dtype=float)
        
        PIP_rjctStartT = np.array([np.nan, np.nan, np.nan],dtype=float)
        PIP_rjctEndT = np.array([np.nan, np.nan, np.nan],dtype=float)

    elif flight == '20150706':
        flFile = '20150706I1_AC.nc'
        dewPtSens = 'TDM.2'
        altSens = 'AltGPS.2'
        wsSens = 'WS.d'

        startT = np.array([11987, 12653, 15824, 16805, 20419, 21384, 22950, 23879],dtype=float)
        endT = np.array([12652, 13534, 16804, 17615, 21349, 22393, 23878, 24875],dtype=float)
        
        PDDstartT = np.array([6988, 10269, 11536, 13570, 14475, 15460, 17760, 19098, 22517, 25204],dtype=float)
        PDDendT = np.array([8340, 11435, 12063, 14172, 15380, 15820, 18897, 20401, 22950, 27063],dtype=float)

        mlBotTemp = np.array([np.nan, 1.9535, 3.7707, 4.9644, 1.9219, 3.5011, 3.0759, 3.6202],dtype=float) 
        mlBotTime = np.array([np.nan, 13099, 16228, 17252, 20851, 21872, 23420, 24383],dtype=float) 
        
        mlTopTemp = np.array([np.nan, 0.8741, 0.0271, 0.2957, 0.0111, 0.3134, 0.1221, 0.0554],dtype=float) 
        mlTopTime = np.array([np.nan, 13054, 16336, 17156, 20917, 21792, 23492, 24281],dtype=float)

        mcsType = np.array(['F','F','T','T','T','T','T','T'],dtype=str)

        sprlZone = np.array(['X', 'X', 'S', 'S', 'S', 'S', 'S', 'S'],dtype=str)

        CIP_intArvThrsh = np.array([1.5e-5, 2.7e-5, 2.7e-5, 3.4e-5, 3.3e-6, 9e-6, 2.5e-5, 2.5e-5],dtype=float)

        PIP_rjctStartT = np.array([np.nan, 12799, np.nan, 17068, np.nan, 21660,
                                  startT[6], startT[7]],dtype=float)
        PIP_rjctEndT = np.array([np.nan, endT[1], np.nan, endT[3], np.nan, endT[5],
                                endT[6], endT[7]],dtype=float)
        
    elif flight == '20150709':
        flFile = '20150709I1_AC.nc'
        dewPtSens = 'TDM.2'
        altSens = 'AltGPS.3'
        wsSens = 'WS.d'

        startT = np.array([8857, 9584, 10643, 11371, 12769, 13665, 14517, 15379, 
                           16585, 17482,19045, 19960, 20891, 21834, 22729, 23546],dtype=float)
        endT = np.array([9582, 10417, 11369, 12160, 13532, 14516, 15377, 16150, 
                         17481, 18226, 19959, 20687, 21833, 22363, 23544, 24074],dtype=float)
        
        PDDstartT = np.array([5975, 7245, 8433, 10404, 12117, 16351, 18316, 18695, 20694, 22451, 24082],dtype=float)
        PDDendT = np.array([7161, 8330, 8847, 10656, 12720, 16602, 18574, 18876, 20885, 22717, 24732],dtype=float)

        mlBotTemp = np.array([0.827, 2.76, 1.584, 2.2527, 2.8054, 4.7628, 3.823, 5.7988, 3.7277, 
                              5.6487, 2.3316, 4.7067, 3.6857, 4.9586, 2.5909, 5.7844],dtype=float)
        mlBotTime = np.array([9166, 10102, 10997, 11769, 13072, 14202, 14867, 15861, 16956, 17932, 
                              19429, 20445, 21312, 22169, 23075, 23942],dtype=float)
        
        mlTopTemp = np.array([0.0104, 0.0034, 0.0423, 0.1119, 0.4757, 3.4426, 1.3057, 2.5534, 0.4398, 
                              2.0643, 0.0656, 1.5315, 0.6095, 2.740, 0.6454, 2.5692],dtype=float)
        mlTopTime = np.array([9185, 10028, 11036, 11711, 13146, 14163, 14903, 15820, 17038, 
                              17866, 19494, 20391, 21391, 22136, 23121, 23882],dtype=float)

        mcsType = np.array(['P','P','P','P','P','P','P','P',
                            'T','T','T','T','T','T','T','T'],dtype=str)

        sprlZone = np.array(['A', 'A', 'A', 'A', 'S', 'S', 'S', 'S', 'S', 
                             'S', 'S', 'S', 'S', 'S', 'S', 'S'],dtype=str)

        CIP_intArvThrsh = np.array([6.9e-5, 5.4e-5 ,6.9e-5, 3.4e-5, 3.4e-5, 2.7e-5, 1e-4, 8.7e-5, 4.3e-5, 
                                    1e-4, 3.4e-5, 3.4e-5, 3.4e-5, 2.7e-5, 3.4e-5, 3.4e-5],dtype=float)
        
        PIP_rjctStartT = np.array([np.nan, 9933, startT[2], 11711, np.nan, 14235, startT[6], 15779, np.nan, 
                                   17845, np.nan, 20372, startT[12], 21986, startT[14], 23678],dtype=float)
        PIP_rjctEndT = np.array([np.nan, endT[1], 10996, endT[3], np.nan, endT[5], 14878, endT[7], np.nan, 
                                   endT[9], np.nan, endT[11], 21120, endT[13], 22818, endT[15]],dtype=float)
        
    else:
        sys.exit('The flight string provided is not valid')

    numSprls = len(startT)
    numPDDs = len(PDDstartT)

    
    ## Define flight-independent parameters
    # Variables used from flight-level data
    latSens = 'LatGPS.1'
    lonSens = 'LonGPS.1'
    tempSens = 'TA.d'
    rhSens = 'HUM_REL.d'
    spSens = 'PSM.2'
    dpSens = 'PQM.2'
    wSens = 'DPJ_WSZ'
    wind_xSens = 'UWX.d'
    wind_ySens = 'UWY.d'
    relWind_xSens = 'URX.d'
    relWind_ySens = 'URY.d'
    wdSens = 'WD.d'
    

    
    ## Create the netCDF file and define global/variable attributes
    rootGrp = Dataset((savePath + flight + '_PECANparams.nc'),'w',format='NETCDF4')
    rootGrp.set_fill_on()

    sprls = rootGrp.createDimension('sprls',numSprls)
    pdds = rootGrp.createDimension('pdds',numPDDs)

    sT = rootGrp.createVariable('startT','f8',('sprls',))
    eT = rootGrp.createVariable('endT','f8',('sprls',))
    PsT = rootGrp.createVariable('PDDstartT','f8',('pdds',))
    PeT = rootGrp.createVariable('PDDendT','f8',('pdds',))
    mlBtmp = rootGrp.createVariable('mlBotTemp','f8',('sprls',))
    mlBtime = rootGrp.createVariable('mlBotTime','f8',('sprls',))
    mlTtmp = rootGrp.createVariable('mlTopTemp','f8',('sprls',))
    mlTtime = rootGrp.createVariable('mlTopTime','f8',('sprls',))
    mcsT = rootGrp.createVariable('mcsType','S1',('sprls',))
    zone = rootGrp.createVariable('sprlZone','S1',('sprls',))
    CiaT = rootGrp.createVariable('CIP_intArvThrsh','f8',('sprls',))
    PRjctST = rootGrp.createVariable('PIP_rjctStartT','f8',('sprls',))
    PRjctET = rootGrp.createVariable('PIP_rjctEndT','f8',('sprls',))

    # Global attributes
    rootGrp.description = 'Flight- and spiral-specific PECAN microphysics parameters'
    rootGrp.flight = flight
    rootGrp.history = 'Created ' + time.asctime(time.gmtime()) + ' UTC'
    rootGrp.FL_rawFile = flFile
    rootGrp.FL_procFile = flight + '_FltLvl_Processed.nc'
    rootGrp.FL_dwptSens = dewPtSens
    rootGrp.FL_altSens = altSens
    rootGrp.FL_latSens = latSens
    rootGrp.FL_lonSens = lonSens
    rootGrp.FL_tempSens = tempSens
    rootGrp.FL_rhSens = rhSens
    rootGrp.FL_spSens = spSens
    rootGrp.FL_dpSens = dpSens
    rootGrp.FL_wSens = wSens
    rootGrp.FL_wind_xSens = wind_xSens
    rootGrp.FL_wind_ySens = wind_ySens
    rootGrp.FL_relWind_xSens = relWind_xSens
    rootGrp.FL_relWind_ySens = relWind_ySens
    rootGrp.FL_wdSens = wdSens
    rootGrp.FL_wsSens = wsSens


    sT.description = 'Time spiral began'
    sT.units = 'Seconds since midnight UTC'
    sT._Fill_Value = np.nan

    eT.description = 'Time spiral ended'
    eT.units = 'Seconds since midnight UTC'
    eT._Fill_Value = np.nan
    
    PsT.description = 'Time PDD began'
    PsT.units = 'Seconds since midnight UTC'
    PsT._Fill_Value = np.nan

    PeT.description = 'Time PDD ended'
    PeT.units = 'Seconds since midnight UTC'
    PeT._Fill_Value = np.nan

    mlBtmp.long_name = 'Melting layer bottom temperature'
    mlBtmp.description = 'Temperature where last evidence of ice was observed in spiral'
    mlBtmp.units = 'deg C'
    mlBtmp._Fill_Value = np.nan

    mlBtime.long_name = 'Melting layer bottom time'
    mlBtime.description = 'Time when last evidence of ice was observed in spiral'
    mlBtime.units = 'Seconds since midnight UTC'
    mlBtime._Fill_Value = np.nan
    
    mlTtmp.long_name = 'Melting layer top temperature'
    mlTtmp.description = 'Temperature where first evidence of liquid water was observed in spiral'
    mlTtmp.units = 'deg C'
    mlTtmp._Fill_Value = np.nan

    mlTtime.long_name = 'Melting layer top time'
    mlTtime.description = 'Time when first evidence of liquid water was observed in spiral'
    mlTtime.units = 'Seconds since midnight UTC'
    mlTtime._Fill_Value = np.nan

    mcsT.description = 'General MCS system type at time of each spiral'
    mcsT.units = 'T=Trailing Stratiform,L=Leading Stratiform,P=Parallel Stratiform,C=Cluster MCS,F=Post-Frontal'
    
    zone.description = 'Location of spiral relative to MCS structure'
    zone.units = 'T = Transition Zone; S = Enhanced Stratiform Region; A = Anvil Region; X = N/A'

    CiaT.description = 'Thresholds of CIP particle interarrival time (in seconds) for each spiral below which particles are to be considered shattered'
    CiaT.units = 'seconds'
    CiaT._Fill_Value = np.nan
    
    PRjctST.description = 'Start time of bad PIP data for each spiral'
    PRjctST.units = 'Seconds since midnight UTC'
    PRjctST._Fill_Value = np.nan
    
    PRjctET.description = 'End time of bad PIP data for each spiral'
    PRjctET.units = 'Seconds since midnight UTC'
    PRjctET._Fill_Value = np.nan

    
    ## Write variables to file
    sT[:] = startT
    eT[:] = endT
    PsT[:] = PDDstartT
    PeT[:] = PDDendT
    mlBtmp[:] = mlBotTemp
    mlBtime[:] = mlBotTime
    mlTtmp[:] = mlTopTemp
    mlTtime[:] = mlTopTime
    mcsT[:] = mcsType
    zone[:] = sprlZone
    CiaT[:] = CIP_intArvThrsh
    
    PRjctST[:] = PIP_rjctStartT
    PRjctET[:] = PIP_rjctEndT

    rootGrp.close()