Blame view

tr_snake.py 7.85 KB
509db604   Paul Nicot   maj tr_position.p...
1
#!/usr/bin/env python
b903d9e8   Paul Nicot   ajout de la trace...
2
# -*- coding: utf-8 -*-
509db604   Paul Nicot   maj tr_position.p...
3

278eaffa   Paul Nicot   ajout de commenta...
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
"""
	Plot maps of the spots tracks, over different parts of the Gulf of St Lawrence.
	Spot positions extract by the script tr_position.py.
		-» map of the entire Gulf of St Lawrence
		-» map of the St Lawrence Estuary
		-» map of the Baie des Chaleurs

	the track are compiled from the selected date to the present, the last 48h are
	ploted with a colorscale to better visualise the spot motion.
		
	Autor: 		Paul Nicot - UQAR-ISMER
	Last updade:	Dec 2015

"""

import matplotlib
matplotlib.use('Agg')
# for running without X server
509db604   Paul Nicot   maj tr_position.p...
22
23
24
25
from mpl_toolkits.basemap import Basemap
from time import gmtime, strftime
import matplotlib.pyplot as plt
import time
c0cb8ae6   Paul Nicot   ajout de la carte...
26
import datetime
509db604   Paul Nicot   maj tr_position.p...
27
28
29
30
import numpy as np
import pandas as pd
from mpl_toolkits.axes_grid1 import make_axes_locatable

278eaffa   Paul Nicot   ajout de commenta...
31
32
# Modifiable parameters
# ================================================
b903d9e8   Paul Nicot   ajout de la trace...
33
snake = 48				# nombre d'heure de suivi du serpent
278eaffa   Paul Nicot   ajout de commenta...
34
start = 20150609				# date de debut de suivi de la saison
b903d9e8   Paul Nicot   ajout de la trace...
35

278eaffa   Paul Nicot   ajout de commenta...
36
37
file_pos = '/share/archives/partage_lasso/spot/spots_pos_flag_clean.csv'
dir_carte = '/share/archives/partage_lasso/spot/carte'
c0cb8ae6   Paul Nicot   ajout de la carte...
38
39
# ================================================

278eaffa   Paul Nicot   ajout de commenta...
40
41
42

# fixed parameters
#-------------------------------------------------
c0cb8ae6   Paul Nicot   ajout de la carte...
43
44
sday  = 84600
shour = 3600
509db604   Paul Nicot   maj tr_position.p...
45
46
47

# colormap
cmap=plt.get_cmap('OrRd')
c0cb8ae6   Paul Nicot   ajout de la carte...
48
49
50
51
52
lab1 = int(snake)
lab2 = int(snake*0.75)
lab3 = int(snake*0.5)
lab4 = int(snake*0.25)
lab5 = 0
509db604   Paul Nicot   maj tr_position.p...
53

278eaffa   Paul Nicot   ajout de commenta...
54
55
# database 
df = pd.read_csv(file_pos,delimiter=';')
509db604   Paul Nicot   maj tr_position.p...
56

278eaffa   Paul Nicot   ajout de commenta...
57
# time paramters
c0cb8ae6   Paul Nicot   ajout de la carte...
58
59
60
61
62
63
64
timestamp = int(time.time())
dt = timestamp - shour * snake
diff = timestamp - dt

strtime = datetime.datetime.fromtimestamp(int(timestamp)).strftime('%Y-%m-%d %H:%M')
nametime  = datetime.datetime.fromtimestamp(int(timestamp)).strftime('%Y%m%d%H%M')

278eaffa   Paul Nicot   ajout de commenta...
65
# date of tracks's beginning
b903d9e8   Paul Nicot   ajout de la trace...
66
67
sT = str(start)
starT = int(time.mktime(datetime.datetime.strptime(sT, "%Y%m%d").timetuple()))
c0cb8ae6   Paul Nicot   ajout de la carte...
68

278eaffa   Paul Nicot   ajout de commenta...
69
70
71
72
73
74

#===================================================
# Gulf of St Lawrence tracks map
#===================================================

# selection of the reliable positions
bc4df3fa   Paul Nicot   maj
75
GSL = df[(df.FLAG == "OCN") & (df.TIME <= timestamp) 
c0cb8ae6   Paul Nicot   ajout de la carte...
76
		& (df.TIME > timestamp - shour*snake)]
b903d9e8   Paul Nicot   ajout de la trace...
77
78
79
GSL_start = df.loc[(df.FLAG == "OCN") & (df.TIME <= timestamp) 
		& (df.TIME > starT)]

509db604   Paul Nicot   maj tr_position.p...
80
81
fig, ax = plt.subplots()

278eaffa   Paul Nicot   ajout de commenta...
82
# map setting
c0cb8ae6   Paul Nicot   ajout de la carte...
83
m = Basemap(projection='stere',lon_0=-64.5,lat_0=48,\
278eaffa   Paul Nicot   ajout de commenta...
84
85
86
			llcrnrlat=45.2,urcrnrlat=51.6,\
			llcrnrlon=-70,urcrnrlon=-56.5,\
			rsphere=6371.2,resolution='h')
509db604   Paul Nicot   maj tr_position.p...
87

278eaffa   Paul Nicot   ajout de commenta...
88
# map background
509db604   Paul Nicot   maj tr_position.p...
89
90
91
92
m.drawmapboundary(fill_color='lightblue') # fill to edge
m.drawcountries()
m.fillcontinents(color='grey',lake_color='lightblue',zorder=0)
m.drawrivers()
c0cb8ae6   Paul Nicot   ajout de la carte...
93
94
95
m.drawcoastlines(linewidth=0.7)
m.drawparallels(np.arange(40,60,2),labels=[1,0,0,0],linewidth=0.0)
m.drawmeridians(np.arange(-80,-40,2),labels=[0,0,0,1],linewidth=0.0)
509db604   Paul Nicot   maj tr_position.p...
96

278eaffa   Paul Nicot   ajout de commenta...
97
# plot track from the start date
b903d9e8   Paul Nicot   ajout de la trace...
98
99
100
101
for spot, group in GSL_start.groupby(['SPOT']):
	lat = group.LAT.values
	lon = group.LON.values
	x,y = m(lon, lat)
278eaffa   Paul Nicot   ajout de commenta...
102
103
104
105
106
107
108
109
110
111
112
113
114
115
	plt.plot(x,y,color='0.35',linewidth=0.7)

# plot colored track of the last positions (snake)
if not GSL.empty:
	for spot, group in GSL.groupby(['SPOT']):
		latitude = group.LAT.values
		longitude = group.LON.values
		t = (group.TIME.values -float(dt)) / diff
		x,y = m(longitude, latitude)
		cax = plt.scatter(x,y,c=t, s=35, cmap=cmap,linewidth=0,vmin=0, vmax=1)

if GSL.empty:
	cax = plt.scatter(-1,-1,c=1, s=35, cmap=cmap,linewidth=0,vmin=0, vmax=1)

509db604   Paul Nicot   maj tr_position.p...
116
117
# colorbar
cbar = fig.colorbar(cax,fraction=0.03, pad=0.03)
c0cb8ae6   Paul Nicot   ajout de la carte...
118
cbar.set_ticks([0,0.25,0.5,0.75,1])
509db604   Paul Nicot   maj tr_position.p...
119
cbar.update_ticks()
c0cb8ae6   Paul Nicot   ajout de la carte...
120
121
cbar.ax.set_yticklabels([lab1,lab2,lab3,lab4,lab5])
cbar.set_label('heures', rotation=270,labelpad=15)
509db604   Paul Nicot   maj tr_position.p...
122
123

# legende
c0cb8ae6   Paul Nicot   ajout de la carte...
124
plt.annotate(strtime + ' HAE', xy=(1, 1),  xycoords='data', size=20,
278eaffa   Paul Nicot   ajout de commenta...
125
126
			xytext=(0.7, 0.9), textcoords='axes fraction',
			horizontalalignment='right', verticalalignment='top')
b903d9e8   Paul Nicot   ajout de la trace...
127
plt.title(u"Dérive de surface (" + str(snake) + u" dernières heures)")
278eaffa   Paul Nicot   ajout de commenta...
128
129
130
# recording
plt.savefig(dir_carte + '/GSL/GSL_' + nametime + '.png',dpi=300)
# only png can be generate
c0cb8ae6   Paul Nicot   ajout de la carte...
131
plt.close()
c0cb8ae6   Paul Nicot   ajout de la carte...
132

278eaffa   Paul Nicot   ajout de commenta...
133
134
135
136
#===================================================
# St Lawrence Estuary
#===================================================
# selection of the reliable positions
b903d9e8   Paul Nicot   ajout de la trace...
137
SLE = df[(df.FLAG == "OCN")
278eaffa   Paul Nicot   ajout de commenta...
138
	   & (df.TIME <= timestamp) & (df.TIME > timestamp - shour*snake)]
b903d9e8   Paul Nicot   ajout de la trace...
139
SLE_start = df.loc[(df.FLAG == "OCN") & (df.TIME <= timestamp)
278eaffa   Paul Nicot   ajout de commenta...
140
			& (df.TIME > starT)]
509db604   Paul Nicot   maj tr_position.p...
141

509db604   Paul Nicot   maj tr_position.p...
142
143
fig, ax = plt.subplots()

278eaffa   Paul Nicot   ajout de commenta...
144
# map setting 
509db604   Paul Nicot   maj tr_position.p...
145
m = Basemap(projection='stere',lon_0=-68.5,lat_0=48,\
278eaffa   Paul Nicot   ajout de commenta...
146
147
148
			llcrnrlat=47.8,urcrnrlat=49.5,\
			llcrnrlon=-70,urcrnrlon=-67,\
			rsphere=6371.2,resolution='h')
509db604   Paul Nicot   maj tr_position.p...
149
150
151
152
153

m.drawmapboundary(fill_color='lightblue') # fill to edge
m.drawcountries()
m.fillcontinents(color='grey',lake_color='lightblue',zorder=0)
m.drawrivers()
c0cb8ae6   Paul Nicot   ajout de la carte...
154
155
156
m.drawcoastlines(linewidth=0.7)
m.drawparallels(np.arange(40,60,0.5),labels=[1,0,0,0],linewidth=0.0)
m.drawmeridians(np.arange(-80,-40,1),labels=[0,0,0,1],linewidth=0.0)
509db604   Paul Nicot   maj tr_position.p...
157

278eaffa   Paul Nicot   ajout de commenta...
158
# plot track from the start date 
b903d9e8   Paul Nicot   ajout de la trace...
159
160
161
162
163
164
for spot, group in SLE_start.groupby(['SPOT']):
	lat = group.LAT.values
	lon = group.LON.values
	x,y = m(lon, lat)
	plt.plot(x,y,color='0.35',linewidth=0.7)

278eaffa   Paul Nicot   ajout de commenta...
165
166
167
168
169
170
171
172
173
174
175
176
# plot track of the last positions (snake)
if not SLE.empty:
	for spot, group in SLE.groupby(['SPOT']):
		latitude = group.LAT.values
		longitude = group.LON.values
		t = (group.TIME.values - float(dt)) / diff
		x,y = m(longitude, latitude)
		cax = plt.scatter(x,y,c=t, s=45, cmap=cmap,linewidth=0,vmin=0, vmax=1)

if SLE.empty:
	cax = plt.scatter(-1,-1,c=1, s=35, cmap=cmap,linewidth=0,vmin=0, vmax=1)	

509db604   Paul Nicot   maj tr_position.p...
177
178
# colorbar
cbar = fig.colorbar(cax,fraction=0.03, pad=0.03)
c0cb8ae6   Paul Nicot   ajout de la carte...
179
cbar.set_ticks([0,0.25,0.5,0.75,1])
509db604   Paul Nicot   maj tr_position.p...
180
cbar.update_ticks()
c0cb8ae6   Paul Nicot   ajout de la carte...
181
182
cbar.ax.set_yticklabels([lab1,lab2,lab3,lab4,lab5])
cbar.set_label('heures', rotation=270,labelpad=15)
509db604   Paul Nicot   maj tr_position.p...
183
184

# legende
c0cb8ae6   Paul Nicot   ajout de la carte...
185
plt.annotate(strtime + ' HAE', xy=(1, 1),  xycoords='data', size=20,
278eaffa   Paul Nicot   ajout de commenta...
186
187
			xytext=(0.9, 0.1), textcoords='axes fraction',
			horizontalalignment='right', verticalalignment='top')
b903d9e8   Paul Nicot   ajout de la trace...
188
plt.title(u"Dérive de surface (" + str(snake) + u" dernières heures)")
278eaffa   Paul Nicot   ajout de commenta...
189
190
# recording
plt.savefig(dir_carte + '/SLE/SLE_' + nametime + '.png',dpi=300)
509db604   Paul Nicot   maj tr_position.p...
191

278eaffa   Paul Nicot   ajout de commenta...
192
plt.close()
b903d9e8   Paul Nicot   ajout de la trace...
193

278eaffa   Paul Nicot   ajout de commenta...
194
195
196
197
198
199
200
201
#===================================================	
# Baie des chaleurs
#===================================================
# selection of the reliable positions
BDC = df[(df.FLAG == "OCN") 
	   & (df.TIME <= timestamp) & (df.TIME > timestamp - shour*snake)]
BDC_start = df.loc[(df.FLAG == "OCN") & (df.TIME <= timestamp)
					& (df.TIME > starT)]
509db604   Paul Nicot   maj tr_position.p...
202
203
fig, ax = plt.subplots()

278eaffa   Paul Nicot   ajout de commenta...
204
# map setting
c0cb8ae6   Paul Nicot   ajout de la carte...
205
m = Basemap(projection='stere',lon_0=-68.5,lat_0=48,\
278eaffa   Paul Nicot   ajout de commenta...
206
207
208
			llcrnrlat=47.5,urcrnrlat=48.9,\
			llcrnrlon=-66.9,urcrnrlon=-63.8,\
			rsphere=6371.2,resolution='h')
509db604   Paul Nicot   maj tr_position.p...
209
210
211
212
213

m.drawmapboundary(fill_color='lightblue') # fill to edge
m.drawcountries()
m.fillcontinents(color='grey',lake_color='lightblue',zorder=0)
m.drawrivers()
c0cb8ae6   Paul Nicot   ajout de la carte...
214
215
216
m.drawcoastlines(linewidth=0.7)
m.drawparallels(np.arange(40,60,0.5),labels=[1,0,0,0],linewidth=0.0)
m.drawmeridians(np.arange(-80,-40,1),labels=[0,0,0,1],linewidth=0.0)
509db604   Paul Nicot   maj tr_position.p...
217

278eaffa   Paul Nicot   ajout de commenta...
218
# plot track from the start date
b903d9e8   Paul Nicot   ajout de la trace...
219
220
221
222
223
224
for spot, group in BDC_start.groupby(['SPOT']):
	lat = group.LAT.values
	lon = group.LON.values
	x,y = m(lon, lat)
	plt.plot(x,y,color='0.35',linewidth=0.7)	

278eaffa   Paul Nicot   ajout de commenta...
225
226
227
228
229
230
231
232
233
234
235
236
# plot track of the last positions (snake)
if not BDC.empty:
	for spot, group in BDC.groupby(['SPOT']):
		latitude = group.LAT.values
		longitude = group.LON.values
		t = (group.TIME.values - float(dt)) / diff
		x,y = m(longitude, latitude)
		cax = plt.scatter(x,y,c=t, s=45, cmap=cmap,linewidth=0,vmin=0, vmax=1)

if BDC.empty:
	cax = plt.scatter(-1,-1,c=1, s=35, cmap=cmap,linewidth=0,vmin=0, vmax=1)	

509db604   Paul Nicot   maj tr_position.p...
237
238
# colorbar
cbar = fig.colorbar(cax,fraction=0.03, pad=0.03)
c0cb8ae6   Paul Nicot   ajout de la carte...
239
cbar.set_ticks([0,0.25,0.5,0.75,1])
509db604   Paul Nicot   maj tr_position.p...
240
cbar.update_ticks()
c0cb8ae6   Paul Nicot   ajout de la carte...
241
242
cbar.ax.set_yticklabels([lab1,lab2,lab3,lab4,lab5])
cbar.set_label('heures', rotation=270,labelpad=15)
509db604   Paul Nicot   maj tr_position.p...
243
244

# legende
c0cb8ae6   Paul Nicot   ajout de la carte...
245
plt.annotate(strtime + ' HAE', xy=(1, 1),  xycoords='data', size=20,
278eaffa   Paul Nicot   ajout de commenta...
246
247
248
			xytext=(0.6, 0.1), textcoords='axes fraction',
			horizontalalignment='right', verticalalignment='top') 
			
b903d9e8   Paul Nicot   ajout de la trace...
249
plt.title(u"Dérive de surface (" + str(snake) + u" dernières heures)")
278eaffa   Paul Nicot   ajout de commenta...
250
251
# recording
plt.savefig(dir_carte + '/BDC/BDC_' + nametime + '.png',dpi=300)
c0cb8ae6   Paul Nicot   ajout de la carte...
252
plt.close()
278eaffa   Paul Nicot   ajout de commenta...