TransPAC_code/hysplit_traceback_single.py at master · dashamstyr/TransPAC_code · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
def traceback(level, station_id, aerofilt_dir):
    import hysplit_tools as tools
    import numpy as np
    import scipy.io
    import os,sys
    import math

    Tman_coord = [39.0,84.0]
    Tman_a = 800 #km semi-major axis
    Tman_b = 300 #km semi-minor axis
    Tman_tilt = 70 #degrees between North and major axis

    Gobi_coord = [43.0,106.0]
    Gobi_a = 1000 #km semi_major axis
    Gobi_b = 500
    Gobi_tilt = 65 #degrees between North and major axis

    G_maxlon = 118
    G_minlon = 95

    G_maxlat = 49
    G_minlat = 37

    T_maxlon = 93
    T_minlon = 75

    T_maxlat = 43
    T_minlat = 35


    #load .mat file

    data_keys = ('station','start_loc','start_date','end_loc','delta_t','delta_d','desert_tag')
    startdir = os.getcwd()
    topdir = aerofilt_dir+'\Level '+level+'\\'+station_id
    os.chdir(topdir)
    output_filename = station_id+'traceback2'

    hysplit_files = os.listdir(os.getcwd())

    data_index = []
    start_date = []
    start_loc = []
    end_loc = []
    delta_t = []
    station = []
    delta_d = []
    desert_tag = []

    print 'Processing %s folder' %station_id

    for h in hysplit_files:

        filetest = h.split('.')

        try:
            filetest[2] == 'mat'

            hysplit_dict = scipy.io.loadmat(h)

            lat = hysplit_dict['lat']
            lon = hysplit_dict['lon']

            data_index = []

            d_total = 0
            for n in range(1,len(lat)):

                [d_temp,theta_temp] = tools.haversine(lat[n-1],lon[n-1],lat[n],lon[n])
                d_total = d_total + d_temp

                if G_minlat < lat[n] < G_maxlat:
                    if G_minlon < lon[n] < G_maxlon:

                        [G_range,G_theta] = tools.haversine(lat[n],lon[n],Gobi_coord[0],Gobi_coord[1])

                        G_rad = tools.ellipserad(Gobi_a,Gobi_b,G_theta,Gobi_tilt)

                        if G_range < G_rad:
                            data_index.append(n)
                            delta_d.append(d_total)
                            desert_tag.append('Gobi')

                if T_minlat < lat[n] < T_maxlat:
                    if T_minlon < lon[n] < T_maxlon:

                        [T_range,T_theta] = tools.haversine(lat[n],lon[n],Tman_coord[0],Tman_coord[1])
                        T_rad = tools.ellipserad(Tman_a,Tman_b,T_theta,Tman_tilt)

                        if T_range < T_rad:
                            data_index.append(n)
                            delta_d.append(d_total)
                            desert_tag.append('Tman')

            for n in range(0,len(data_index)):
                start_date.append((hysplit_dict['year'][0],hysplit_dict['month'][0],hysplit_dict['day'][0],\
                                 hysplit_dict['hour'][0]))
                start_loc.append((lat[0],lon[0]))
                end_loc.append((lat[data_index[n]],lon[data_index[n]]))
                delta_t.append(hysplit_dict['delta_t'][data_index[n]])
                station.append(station_id)

        except IndexError:
            pass

    output_data = (station,start_loc,start_date,end_loc,delta_t,delta_d,desert_tag)

    output_dict = dict(zip(data_keys,output_data))

    scipy.io.savemat(output_filename,output_dict)


    os.chdir(startdir)