source: branches/2017/dev_v3.20_2017_transport_max/GRAPHICS/ISPOL_profiles.pro

Last change on this file was 6, checked in by vancop, 8 years ago

initial import of v3.20 /Users/ioulianikolskaia/Boulot/CODES/LIM1D/ARCHIVE/TMP/LIM1D_v3.20/
File size: 26.1 KB
Line 
1;
2;------------------------------------------------------------------------------
3;
4; This script makes plots of biological simulations
5;
6; chla stock
7; NO3 stock
8; PO4 stock
9;
10; (c) Martin Vancoppenolle, UCL-ASTR, june 2007
11; reviewed may 2008
12; reviewed feb 2010, for biological model
13; reviewed aug 2014, helsinki
14;
15;------------------------------------------------------------------------------
16;
17miss_val = -9999.
18i_auto = 1             ; 1=use exp_id.dat; 2=use manual mode with sensitivity runs
19plot_type = 'DEFAULT'  ;
20;
21;==============================================================================
22; General Options
23;==============================================================================
24;
25
26IF ( i_auto EQ 1 ) THEN BEGIN ; automatic mode
27
28   nruns = 1
29   control_run = 'E001' & indir = '' & outdir = '' & dt = 86400 & c_bio_model = 'BFMSI' & site = 'TOURNAI'
30   
31   OPENR, 11, 'exp_id.dat'
32   READF, 11, control_run & READF, 11, indir & READF, 11, outdir
33   READF, 11, dt & READF, 11, c_bio_model & READF, 11, site
34   CLOSE, 11
35   
36   PRINT, control_run
37   PRINT, indir
38   PRINT, outdir
39   PRINT, dt
40   PRINT, c_bio_model
41   PRINT, site
42
43ENDIF
44
45;-------------
46; Manual mode
47;-------------
48IF ( i_auto EQ 0 ) THEN BEGIN
49
50   indir_base = '/Users/ioulianikolskaia/Boulot/CODES/BFMLIM_2013/BFM-LIM-2506_2014/LIM/RUN/'
51   outdir = '/Users/ioulianikolskaia/Boulot/SCIENCE/PLOT_SCRIPTS/LIM1D_BIO/IDL/plots/'
52   dt   = 3600.        ; time step
53
54   control_run = 'LIM_KR_ML'
55   sensiti_run = [ ]
56;  sensiti_run =  [ 'LIM_KR_SL', 'LIM_KR_BA' ]
57   dummy = SIZE(sensiti_run)
58   nruns = dummy(1) + 1 ; sensitivity runs + 1 control run
59   print, 'nruns ', nruns
60
61   indir = STRARR(nruns)
62   indir(0) = indir_base+control_run+'/'
63   FOR i = 1, nruns-1 DO BEGIN
64      indir(i) = indir_base+sensiti_run(i-1)+'/'
65   ENDFOR
66
67ENDIF
68
69;
70;==============================================================================
71; Graphic Options
72;==============================================================================
73;
74
75;--------------------------------
76;--- plot_type dependent options
77;--------------------------------
78
79IF ( plot_type EQ 'DEFAULT' ) THEN BEGIN
80
81   colors      = [  0, 200 ,   100   ,   240   ,   240   ,   120   ,   120    ]     
82   linestyle   = [  0,0    ,    0    ,    1    ,     2   ,     2   ,     1    ]
83   thick       = [  5, 5.  ,    5.0  ,    2.0  ,   2.0   ,   2.0   ,    2.0   ]
84   stamp_out   = control_run+'_profiles'
85   cs          = 1.8 ; charsize
86   colorkey    = 'rd'
87   ct          = 0
88   !X.MARGIN   = [8,3]
89   !Y.MARGIN   = [3,3]
90   x_size      = 4.
91   y_size      = 6.
92   add_obs     = 'YES'
93
94ENDIF
95
96file_name = stamp_out+'.ps'
97
98numplot_x = 5       ; number of horizontal plots
99numplot_y = 2       ; number of vertical plots
100
101device = 'PS'       ; 'PS' or 'X'
102
103;--------------------------------
104;--- define output plots
105;--------------------------------
106figuresize_x = numplot_x * x_size ; figure size on x direction
107figuresize_y = numplot_y * y_size ; figure size on y direction
108set_plot, device
109IF ( device EQ 'PS' ) THEN BEGIN
110   device, /COLOR, /LANDSCAPE, filename=outdir+file_name, $
111           XSIZE=figuresize_x,YSIZE=figuresize_y,FONT_SIZE=9.0
112   loadct, ct
113ENDIF
114
115IF ( device EQ 'X' ) THEN BEGIN
116   xsize = 1200
117   ysize = 800
118   init_graphics_x, xsize, ysize, colorkey
119ENDIF
120
121;
122;==============================================================================
123; EXTRACT DATA
124;==============================================================================
125;
126;open_fields_v6, indir(0), control_run, c_bio_model, $
127;                numt, doy, ts_d, ts_m,                        $
128;                h_i, h_s, z_ip, z_ib, s_i, t_i, e_i, PAR, Ra, $
129;                dhib, dhisu, dhisi,                           $
130;                DAFb, NO3b, PO4b, DSib, chla,                 $
131;                DAFt, NO3t, PO4t, DSit,                       $
132;                FDAbd, FDAbpos, FDAsi,                        $
133;                FNO3, FNO3bpos, FNO3si,                       $
134;                syn, lys, rem,                                $
135;                lim_lig, lim_no3, lim_po4, lim_dsi, lim_tem, lim_sal
136               
137 open_fields_v7, indir(0), control_run, c_bio_model, $
138                 numt, doy, ts_d, ts_m,                        $
139                 h_i, h_s, z_ip, z_ib, s_i, t_i, e_i, PAR, Ra, $
140                 dhib, dhisu, dhisi,                           $
141                 DAFb, NO3b, PO4b, DSib, chla, eoCb,           $
142                 Argb, Argbub, $
143                 Oxyb, Oxybub, $
144                 DICb, Alkb, CO2b, CO2bub,                     $
145                 CO2aq, HCO3m, CO32m, pH, pCO2, Ikab,          $
146                 dFeb, aFeb, eFeb,                             $
147                 DAFt, NO3t, PO4t, DSit, eoCt, Argt,           $
148                 DICt, Alkt, CO2t, Ikat,                       $
149                 dFet, aFet, eFet,                             $
150                 FDAbd, FDAbpos, FDAsi,                        $
151                 FNO3, FNO3bpos, FNO3si,                       $
152                 FCO2_atm, FCO2_bub,                           $
153                 syn, lys, rem,                                $
154                 lim_lig, lim_no3, lim_po4, lim_dsi, lim_tem, lim_sal
155
156zsize = SIZE(DAFb)
157nlay = zsize(1)   ; number of layers
158nts  = zsize(2)   ; number of time steps
159
160;
161;------------------------------------------------------------------------------
162; OBS ...
163;------------------------------------------------------------------------------
164;
165
166; observed normalized profiles from Lannuzel et al. (2007)
167
168N_obs = 7 ; station numbers
169N_dpt = 6 ; number of depths
170
171; DOY
172doy_obs = [  334,    339,    344,    349,   354,   360,  365 ]
173
174; Ice thickness / snow depth (ULB SITE)
175hi_obs = [ 0.90 ,  0.92 ,  0.95 ,  0.87 , 0.94 , 0.94 , 0.83 ]
176hs_obs = [ 0.095,  0.130,  0.25 ,  0.15 , 0.110, 0.088, 0.06 ]
177
178; Depths
179depth_obs = FLTARR(N_obs, N_dpt)
180depth_obs[0,*] = [ 6, 12, 43, 63, 81, 87 ]  / 100. / hi_obs(0)
181depth_obs[1,*] = [ 6, 12, 43, 63, 81, 87 ]  / 100. / hi_obs(1)
182depth_obs[2,*] = [ 6, 12, 43, 63, 81, 87 ]  / 100. / hi_obs(2)
183depth_obs[3,*] = [4.5,12, 43, 63, 77, 83 ]  / 100. / hi_obs(3)
184depth_obs[4,*] = [ 7, 13, 43, 63, 80, 86 ]  / 100. / hi_obs(4)
185depth_obs[5,*] = [ 6, 12, 43, 63, 78, 87 ]  / 100. / hi_obs(5)
186depth_obs[6,*] = [ 6,12.5,43, 62, 77, 83 ]  / 100. / hi_obs(6)
187
188; Chlorophyll
189chla_obs = FLTARR(N_obs, N_dpt)
190chla_obs[0,*] = [ 0.35, 0.41, 0.18, 0.08, 3.98, 23.57] ; same
191chla_obs[1,*] = [ 0.44, 0.44, 0.45, 0.16, 3.44, 26.47] ; same
192chla_obs[2,*] = [ 0.97, 2.22, 0.40, 0.25, 2.56, 21.48] ; not the same
193chla_obs[3,*] = [ 0.69, 1.23, 0.27, 0.30, 2.19, 24.23] ; same
194chla_obs[4,*] = [ 0.78, 1.06, 0.20, 0.65, 4.20, 28.41] ; same as bug fix 1.0
195chla_obs[5,*] = [ 0.70, 0.66, 0.69, 0.29, 3.30, 16.33] ; same
196chla_obs[6,*] = [ 1.69, 0.97, 0.59, 0.95, 3.62, 24.77] ; same
197
198; Silicates
199dsi_obs = FLTARR(N_obs, N_dpt)
200dsi_obs[0,*] = [ 6.1, 1.9, 2.4, 1.5, 2.0,11.2 ]
201dsi_obs[1,*] = [ 0.3, 0.7, 0.7, 0.7, 2.8, 8.1 ]
202dsi_obs[2,*] = [ 2.2, 1.8, 1.5, 1.3, 0.8, 8.9 ]
203dsi_obs[3,*] = [ 0.9, 1.9, 1.2, 0.6, 1.8, 9.9 ]
204dsi_obs[4,*] = [ 0.5, 0.7, 1.6, 0.4, 0.7, 7.8 ]
205dsi_obs[5,*] = [ 0.6, 1.2, 1.7, 1.2, 1.5, 7.6 ]
206dsi_obs[6,*] = [ 0.9, 1.5, 1.4, 2.0, 2.6, miss_val ]
207
208; NOX
209nox_obs = FLTARR(N_obs, N_dpt)
210nox_obs[0,*] = [ 1.2, 0.5, 0.22, 0.3, 0.4, 7.9 ]
211nox_obs[1,*] = [ 0.4, 0.5,  0.5, 0.5, 0.8, 7.4 ]
212nox_obs[2,*] = [ 1.4, 1.3,  0.4, 0.9, 0.5, 7.3 ]
213nox_obs[3,*] = [ 1.4, 1.9,  0.7, 0.3, 1.9, 11.7]
214nox_obs[4,*] = [ 0.8, 1.4,  1.7, 0.7, 0.24, 6.9]
215nox_obs[5,*] = [ 0.3, 0.7,  0.0, 0.0, 0.0, 3.0 ]
216nox_obs[6,*] = [ 0.0, 0.6,  0.3, 0.4, 1.2, miss_val ]
217
218; PO4
219po4_obs = FLTARR(N_obs, N_dpt)
220po4_obs[0,*] = [ 1.0, 0.3, 0.2, 1.3, 0.3, 3.0 ]
221po4_obs[1,*] = [ 0.0, 0.1, 0.15,0.2, 0.2, 2.8 ]
222po4_obs[2,*] = [ 0.1, 0.0, 0.0,0.05, 0.1, 3.2 ]
223po4_obs[3,*] = [ 0.4, 0.5, 0.4, 0.3, 0.6, 3.22]
224po4_obs[4,*] = [ 0.4, 0.6, 0.5, 0.4, 0.4, 2.2 ]
225po4_obs[5,*] = [ 0.1, 0.1, 0.0, 0.0, 0.1, 1.5 ]
226po4_obs[6,*] = [ 0.0, 0.1, 0.1, 0.0, 0.1, miss_val ]
227
228; physical observations
229; T
230t_i_obs = FLTARR(N_obs, N_dpt)
231t_i_obs[0,*] = [ -3.1, -2.4, -2.3, -2.1, -2, -1.9 ]
232t_i_obs[1,*] = [ -1.7, -2.1, -1.9, -1.6, -1.8, -1.8 ]
233t_i_obs[2,*] = [ -1.2, -1.4, -1.9, -1.8, -1.8, -1.8 ]
234t_i_obs[3,*] = [ -1.5, -1.4, -1.7, -1.7, -1.8, -1.9 ]
235t_i_obs[4,*] = [ -0.4, -1.3, -1.4, -1.5, -1.6, -1.6 ]
236t_i_obs[5,*] = [ -0.2, -0.7, -1.2, -1.2, -1.2, -1.4 ]
237t_i_obs[6,*] = [ -1.1, -1.2, -1.4, -1.4, -1.6, -1.8 ]
238
239
240
241; S-physical core
242s_i_obs = FLTARR(N_obs, N_dpt)
243
244s_i_obs[0,*] = [ 8.74000   ,  7.78000    , 6.12000    , 4.60000   ,  5.80000  ,   9.43000 ]
245s_i_obs[1,*] = [ 5.91100   ,  5.58700    , 5.32000    , 3.39600   ,  5.22900  ,   8.50800 ]
246s_i_obs[2,*] = [ 3.62000   ,  4.34000    , 5.53000    , 5.01000   ,  3.99000  ,   5.16000 ]
247s_i_obs[3,*] = [ 8.36000   ,  5.21000    , 3.89000    , 3.39600   ,  4.69000  ,   8.30377 ]
248s_i_obs[4,*] = [ 2.58000   ,  5.97000    , 4.47000    , 3.25000   ,  3.35000  ,   4.65000 ]
249s_i_obs[5,*] = [ 0.67000   ,  1.99000    , 4.99000    , 3.71000   ,  3.54000  ,   4.62000 ]
250s_i_obs[6,*] = [ 4.78000   ,  4.00000    , 3.93000    , 2.82000   ,  3.95000  ,   6.80384 ]
251
252e_i_obs = - 0.054*s_i_obs / t_i_obs
253
254; brine concentrations
255iaddr = where ( dsi_obs NE miss_val )
256dsi_br_obs = FLTARR(N_obs, N_dpt)
257dsi_br_obs(iaddr) = dsi_obs(iaddr) / e_i_obs(iaddr)
258
259iaddr = where ( po4_obs NE miss_val )
260po4_br_obs = FLTARR(N_obs, N_dpt)
261po4_br_obs(iaddr) = po4_obs(iaddr) / e_i_obs(iaddr)
262
263iaddr = where ( nox_obs NE miss_val )
264nox_br_obs = FLTARR(N_obs, N_dpt)
265nox_br_obs(iaddr) = nox_obs(iaddr) / e_i_obs(iaddr)
266
267;*** Mean observed Profiles
268t_i_obs_mean = MEAN( t_i_obs, dimension = 1 )
269t_i_obs_std = STDDEV( t_i_obs, dimension = 1 )
270
271s_i_obs_mean = MEAN( s_i_obs, dimension = 1 )
272s_i_obs_std = STDDEV( s_i_obs, dimension = 1 )
273
274e_i_obs_mean = MEAN( e_i_obs, dimension = 1 )
275e_i_obs_std = STDDEV( e_i_obs, dimension = 1 )
276
277depth_obs_mean = MEAN( depth_obs, DIMENSION = 1 )
278
279chla_obs_mean  = MEAN( chla_obs , DIMENSION = 1 )
280chla_obs_std = STDDEV( chla_obs, dimension = 1 )
281
282dsi_obs_mean = FLTARR(N_dpt) & dsi_obs_std = FLTARR(N_dpt)
283nox_obs_mean = FLTARR(N_dpt) & nox_obs_std = FLTARR(N_dpt)
284po4_obs_mean = FLTARR(N_dpt) & po4_obs_std = FLTARR(N_dpt)
285
286FOR i = 0, N_dpt - 1 DO BEGIN
287   i_addr = WHERE( dsi_obs(*,i) NE miss_val )
288   dsi_obs_mean(i) =  MEAN( dsi_obs(i_addr,i) )
289   dsi_obs_std(i)  = STDDEV( dsi_obs(i_addr,i) )
290   nox_obs_mean(i) =  MEAN( nox_obs(i_addr,i) )
291   nox_obs_std(i)  = STDDEV( nox_obs(i_addr,i) )
292   po4_obs_mean(i) =  MEAN( po4_obs(i_addr,i) )
293   po4_obs_std(i)  = STDDEV( po4_obs(i_addr,i) )
294ENDFOR
295
296;*** Corresponding model mean profiles
297t_i_mod_mean = FLTARR(nlay)  & t_i_mod_mean(*) = 0.
298t_i_mod_std  = FLTARR(nlay)  & t_i_mod_std (*) = 0.
299s_i_mod_mean = FLTARR(nlay)  & s_i_mod_mean(*) = 0.
300s_i_mod_std  = FLTARR(nlay)  & s_i_mod_std (*) = 0.
301e_i_mod_mean = FLTARR(nlay)  & e_i_mod_mean(*) = 0.
302e_i_mod_std  = FLTARR(nlay)  & e_i_mod_std (*) = 0.
303chla_mod_mean = FLTARR(nlay) & chla_mod_mean(*) = 0.
304chla_mod_std  = FLTARR(nlay) & chla_mod_std (*) = 0.
305nox_mod_mean = FLTARR(nlay)  & nox_mod_mean(*) = 0.
306nox_mod_std  = FLTARR(nlay)  & nox_mod_std (*) = 0.
307po4_mod_mean = FLTARR(nlay)  & po4_mod_mean(*) = 0.
308po4_mod_std  = FLTARR(nlay)  & po4_mod_std (*) = 0.
309dsi_mod_mean = FLTARR(nlay)  & dsi_mod_mean(*) = 0.
310dsi_mod_std  = FLTARR(nlay)  & dsi_mod_std (*) = 0.
311FOR i = 0, nlay - 1 DO BEGIN
312
313   FOR i_obs = 0, N_obs - 1 DO BEGIN
314      i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
315      t_i_mod_mean(i)  = t_i_mod_mean(i) + t_i(i,i_mod) / FLOAT(N_obs)
316      s_i_mod_mean(i)  = s_i_mod_mean(i) + s_i(i,i_mod) / FLOAT(N_obs)
317      e_i_mod_mean(i)  = e_i_mod_mean(i) + e_i(i,i_mod) / FLOAT(N_obs)
318      chla_mod_mean(i) = chla_mod_mean(i) + chla(i,i_mod) / FLOAT(N_obs)
319      dsi_mod_mean(i)  = dsi_mod_mean(i) + dsib(i,i_mod) / FLOAT(N_obs)
320      nox_mod_mean(i)  = nox_mod_mean(i) + no3b(i,i_mod) / FLOAT(N_obs)
321      po4_mod_mean(i)  = po4_mod_mean(i) + po4b(i,i_mod) / FLOAT(N_obs)
322   ENDFOR
323
324   zsum   = 0
325   zsum_s = 0
326   zsum_e = 0
327   zsum_c = 0
328   zsum_d = 0
329   zsum_n = 0
330   zsum_p = 0
331   FOR i_obs = 0, N_obs - 1 DO BEGIN
332      i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
333      zsum   = zsum + ( t_i(i,i_mod) - t_i_mod_mean(i) ) ^2 / FLOAT(N_obs)
334      zsum_s = zsum_s + ( s_i(i,i_mod) - s_i_mod_mean(i) ) ^2 / FLOAT(N_obs)
335      zsum_e = zsum_e + ( e_i(i,i_mod) - e_i_mod_mean(i) ) ^2 / FLOAT(N_obs)
336      zsum_c = zsum_c + ( chla(i,i_mod) - chla_mod_mean(i) ) ^2 / FLOAT(N_obs)
337      zsum_d = zsum_d + ( dsib(i,i_mod) - dsi_mod_mean(i) ) ^2 / FLOAT(N_obs)
338      zsum_n = zsum_n + ( no3b(i,i_mod) - nox_mod_mean(i) ) ^2 / FLOAT(N_obs)
339      zsum_p = zsum_p + ( po4b(i,i_mod) - po4_mod_mean(i) ) ^2 / FLOAT(N_obs)
340   ENDFOR
341   t_i_mod_std(i)  = SQRT( zsum   )
342   s_i_mod_std(i)  = SQRT( zsum_s )
343   e_i_mod_std(i)  = SQRT( zsum_e )
344   chla_mod_std(i) = SQRT( zsum_c )
345   nox_mod_std(i)  = SQRT( zsum_n )
346   po4_mod_std(i)  = SQRT( zsum_p )
347   dsi_mod_std(i)  = SQRT( zsum_d )
348
349ENDFOR
350
351;-------
352
353
354;
355;==============================================================================
356; DIAGNOSTICS
357;==============================================================================
358;
359phi=findgen(32)*(!PI*2/32.)
360phi = [ phi, phi(0) ]
361usersym, cos(phi), sin(phi), /fill
362
363;
364;==============================================================================
365; PLOTS
366;==============================================================================
367;
368;
369;==============================================================================
370; Time series of profiles
371;==============================================================================
372;
373;---------
374;--- T ---
375;---------
376!P.MULTI=[0,numplot_x, numplot_y]
377FOR i_obs = 0, N_obs - 1 DO BEGIN
378
379   ; prepare plot
380   LOADCT, 0
381   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
382   PLOT, [ -3, 0. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'deg C', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2
383   XYOUTS, -1, -0.1, ztitle
384
385   ; add obs
386   OPLOT, t_i_obs(i_obs,*), - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
387
388   ; add model
389   i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
390   zz = FLTARR(nlay)
391   FOR i = 0, nlay - 1 DO BEGIN
392      zz(i) = z_ip(i,i_mod) / h_i(i_mod)
393   ENDFOR
394   LOADCT, 3
395   OPLOT, t_i(*,i_mod), -zz, THICK = 3, COLOR = 150
396
397ENDFOR
398
399
400;---------
401;--- S   
402;---------
403!P.MULTI=[0,numplot_x, numplot_y]
404FOR i_obs = 0, N_obs - 1 DO BEGIN
405
406   ; prepare plot
407   LOADCT, 0
408   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
409   PLOT, [0., 15. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'g/kg', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2
410   XYOUTS, 7., -0.1, ztitle
411
412   ; add obs
413   OPLOT, s_i_obs(i_obs,*), - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
414
415   ; add model
416   i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
417   zz = FLTARR(nlay)
418   FOR i = 0, nlay - 1 DO BEGIN
419      zz(i) = z_ip(i,i_mod) / h_i(i_mod)
420   ENDFOR
421   LOADCT, 3
422   OPLOT, s_i(*,i_mod), -zz, THICK = 3, COLOR = 150
423
424ENDFOR
425
426
427;---------
428;--- e   
429;---------
430!P.MULTI=[0,numplot_x, numplot_y]
431FOR i_obs = 0, N_obs - 1 DO BEGIN
432
433   ; prepare plot
434   LOADCT, 0
435   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
436   PLOT, [0., 50. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = '%', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2
437   XYOUTS, 15., -0.1, ztitle
438
439   ; add obs
440   OPLOT, e_i_obs(i_obs,*)*100., - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
441
442   ; add model
443   i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
444   zz = FLTARR(nlay)
445   FOR i = 0, nlay - 1 DO BEGIN
446      zz(i) = z_ip(i,i_mod) / h_i(i_mod)
447   ENDFOR
448   LOADCT, 3
449   OPLOT, e_i(*,i_mod), -zz, THICK = 3, COLOR = 150
450
451ENDFOR
452
453
454;--------------
455;--- chla
456;--------------
457!P.MULTI=[0,numplot_x, numplot_y]
458FOR i_obs = 0, N_obs - 1 DO BEGIN
459
460   ; prepare plot
461   LOADCT, 0
462   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
463   PLOT, [ 0.01, 40. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mg chla / m3', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2
464   XYOUTS, 27., -0.1, ztitle
465
466   ; add obs
467   OPLOT, chla_obs(i_obs,*), - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
468
469   ; add model
470   i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
471   zz = FLTARR(nlay)
472   FOR i = 0, nlay - 1 DO BEGIN
473      zz(i) = z_ib(i,i_mod) / h_i(i_mod)
474   ENDFOR
475   LOADCT, 9
476   OPLOT, chla(*,i_mod), -zz, THICK = 3, COLOR = 150
477
478ENDFOR
479
480
481;------------
482;--- nox 
483;------------
484!P.MULTI=[0,numplot_x, numplot_y]
485
486FOR i_obs = 0, N_obs - 1 DO BEGIN
487
488   ; prepare plot
489   LOADCT, 0
490   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
491   PLOT, [ 0.00, 12. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mmmol N / m3', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2
492   XYOUTS, 2.2, -0.3, ztitle
493
494   ; add obs
495   OPLOT, nox_obs(i_obs,*), - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
496
497   ; add model
498   i_mod = WHERE(doy+1 EQ doy_obs(i_obs))
499   zz = FLTARR(nlay)
500   FOR i = 0, nlay - 1 DO BEGIN
501      zz(i) = z_ib(i,i_mod) / h_i(i_mod)
502   ENDFOR
503   LOADCT, 9
504   OPLOT, no3b(*,i_mod), -zz, THICK = 3, COLOR = 150
505
506   PRINT, ' i_mod : ', i_mod
507   PRINT, ' NO3b  : ', no3b(*,i_mod)
508
509ENDFOR
510
511
512;-------------
513;--- po4 
514;-------------
515!P.MULTI=[0,numplot_x, numplot_y]
516
517FOR i_obs = 0, N_obs - 1 DO BEGIN
518
519   ; prepare plot
520   LOADCT, 0
521   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
522   PLOT, [ 0.00, 4.0 ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mmmol P / m3', XTICKFORMAT='(I2)', XMINOR = 2, YMINOR = 2
523   XYOUTS, 0.37, -0.3, ztitle
524
525   ; add obs
526   OPLOT, po4_obs(i_obs,*), - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
527
528   ; add model
529   i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
530   zz = FLTARR(nlay)
531   FOR i = 0, nlay - 1 DO BEGIN
532      zz(i) = z_ib(i,i_mod) / h_i(i_mod)
533   ENDFOR
534   LOADCT, 9
535   OPLOT, po4b(*,i_mod), -zz, THICK = 3, COLOR = 150
536
537ENDFOR
538
539
540
541;--------------
542;--- DSi 
543;--------------
544!P.MULTI=[0,numplot_x, numplot_y]
545
546FOR i_obs = 0, N_obs - 1 DO BEGIN
547
548   ; prepare plot
549   LOADCT, 0
550   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
551   PLOT, [ 0.00, 12. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mmmol Si / m3', XTICKFORMAT='(F4.1)', XMINOR = 2, YMINOR = 2
552   XYOUTS, 0.37, -0.3, ztitle
553
554   ; add obs
555   OPLOT, dsi_obs(i_obs,*), - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
556
557   ; add model
558   i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
559   zz = FLTARR(nlay)
560   FOR i = 0, nlay - 1 DO BEGIN
561      zz(i) = z_ib(i,i_mod) / h_i(i_mod)
562   ENDFOR
563   LOADCT, 9
564   OPLOT, dsib(*,i_mod), -zz, THICK = 3, COLOR = 150
565
566ENDFOR
567
568;
569;==============================================================================
570; Time series of profiles, brine concentrations
571;==============================================================================
572;
573
574;--- nox_br
575!P.MULTI=[0,numplot_x, numplot_y]
576
577FOR i_obs = 0, N_obs - 1 DO BEGIN
578
579   ; prepare plot
580   LOADCT, 0
581   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
582   PLOT, [ 0., 30. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mmmol N / m3', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2
583   XYOUTS, 20.0, -0.15, ztitle
584
585   ; add obs
586   OPLOT, nox_br_obs(i_obs,*), - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
587
588   ; add model
589   i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
590   zz = FLTARR(nlay)
591   FOR i = 0, nlay - 1 DO BEGIN
592      zz(i) = z_ib(i,i_mod) / h_i(i_mod)
593   ENDFOR
594   LOADCT, 9
595   OPLOT, no3b(*,i_mod)/e_i(*,i_mod)*100., -zz, THICK = 3, COLOR = 150
596   LOADCT, 1
597
598   OPLOT, [ 1.6, 1.6 ], [ -1., 0.], linestyle = 1, color = 150 , thick = 3
599
600ENDFOR
601
602;--- po4_br
603!P.MULTI=[0,numplot_x, numplot_y]
604
605FOR i_obs = 0, N_obs - 1 DO BEGIN
606
607   ; prepare plot
608   LOADCT, 0
609   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
610   PLOT, [ 0.0, 15. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, $
611         SUBTITLE = 'mmmol P / m3', XTICKFORMAT='(I3)', XMINOR = 2, XTICKS = 3 , YMINOR = 2
612   XYOUTS, 10., -0.15, ztitle
613
614   ; add obs
615   OPLOT, po4_br_obs(i_obs,*), - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
616
617   ; add model
618   i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
619   zz = FLTARR(nlay)
620   FOR i = 0, nlay - 1 DO BEGIN
621      zz(i) = z_ib(i,i_mod) / h_i(i_mod)
622   ENDFOR
623   LOADCT, 9
624   OPLOT, po4b(*,i_mod)/e_i(*,i_mod)*100., -zz, THICK = 3, COLOR = 150
625   LOADCT, 1
626
627   OPLOT, [ 0.24, 0.24 ], [ -1., 0.], linestyle = 1, color = 150 , thick = 3
628
629ENDFOR
630
631;--- dsi_br
632!P.MULTI=[0,numplot_x, numplot_y]
633
634FOR i_obs = 0, N_obs - 1 DO BEGIN
635
636   ; prepare plot
637   LOADCT, 0
638   ztitle = STRCOMPRESS(STRING(doy_obs(i_obs)), /REMOVE_ALL)
639   PLOT, [ 0., 30. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mmmol Si / m3', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2
640   XYOUTS, 20.0, -0.15, ztitle
641
642   ; add obs
643   OPLOT, dsi_br_obs(i_obs,*), - depth_obs(i_obs,*), PSYM = 1, THICK = 3, SYMSIZE = 0.8
644
645   ; add model
646   i_mod = WHERE(doy EQ doy_obs(i_obs)-1)
647   zz = FLTARR(nlay)
648   FOR i = 0, nlay - 1 DO BEGIN
649      zz(i) = z_ib(i,i_mod) / h_i(i_mod)
650   ENDFOR
651   LOADCT, 9
652   OPLOT, dsib(*,i_mod)/e_i(*,i_mod)*100., -zz, THICK = 3, COLOR = 150
653   LOADCT, 1
654
655   OPLOT, [ 3.90, 3.90 ], [ -1., 0.], linestyle = 1, color = 150 , thick = 3
656
657ENDFOR
658
659;
660;==============================================================================
661; Mean profiles
662;==============================================================================
663;
664!P.MULTI=[0,numplot_x, numplot_y]
665
666;--- T ---
667
668; prepare plot
669LOADCT, 0
670ztitle = "Mean T"
671PLOT, [ -3, 0. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'deg C', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2, TITLE = ztitle
672
673; add obs
674OPLOT, t_i_obs_mean, - depth_obs_mean, PSYM = 8, THICK = 3, SYMSIZE = 0.8
675FOR i = 0, N_dpt - 1 DO BEGIN
676   OPLOT, [ t_i_obs_mean(i) - t_i_obs_std(i), t_i_obs_mean(i) + t_i_obs_std(i) ], $
677          [ -depth_obs_mean(i), -depth_obs_mean(i) ]
678ENDFOR
679LOADCT, 3
680OPLOT, t_i_mod_mean, -zz, THICK = 3, COLOR = 150
681LOADCT, 0
682OPLOT, t_i_mod_mean+t_i_mod_std, -zz, THICK = 1, COLOR = 150
683OPLOT, t_i_mod_mean-t_i_mod_std, -zz, THICK = 1, COLOR = 150
684
685;--- S ---
686
687; prepare plot
688LOADCT, 0
689ztitle = "Mean S"
690   PLOT, [0., 10. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'g/kg', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2, TITLE = ztitle
691
692; add obs
693OPLOT, s_i_obs_mean, - depth_obs_mean, PSYM = 8, THICK = 3, SYMSIZE = 0.8
694FOR i = 0, N_dpt - 1 DO BEGIN
695   OPLOT, [ s_i_obs_mean(i) - s_i_obs_std(i), s_i_obs_mean(i) + s_i_obs_std(i) ], $
696          [ -depth_obs_mean(i), -depth_obs_mean(i) ]
697ENDFOR
698
699; add model
700LOADCT, 3
701OPLOT, s_i_mod_mean, -zz, THICK = 3, COLOR = 150
702LOADCT, 0
703OPLOT, s_i_mod_mean+s_i_mod_std, -zz, THICK = 1, COLOR = 150
704OPLOT, s_i_mod_mean-s_i_mod_std, -zz, THICK = 1, COLOR = 150
705
706;--- e ---
707
708; prepare plot
709LOADCT, 0
710ztitle = "Mean e"
711   PLOT, [0., 50. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = '%', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2, TITLE = ztitle
712
713; add obs
714OPLOT, e_i_obs_mean*100., - depth_obs_mean, PSYM = 8, THICK = 3, SYMSIZE = 0.8
715FOR i = 0, N_dpt - 1 DO BEGIN
716   OPLOT, [ e_i_obs_mean(i) - e_i_obs_std(i), e_i_obs_mean(i) + e_i_obs_std(i) ] * 100., $
717          [ -depth_obs_mean(i), -depth_obs_mean(i) ]
718ENDFOR
719
720; add model
721LOADCT, 3
722OPLOT, e_i_mod_mean, -zz, THICK = 3, COLOR = 150
723LOADCT, 0
724OPLOT, e_i_mod_mean+e_i_mod_std, -zz, THICK = 1, COLOR = 150
725OPLOT, e_i_mod_mean-e_i_mod_std, -zz, THICK = 1, COLOR = 150
726
727;--- chl-a ---
728
729; prepare plot
730LOADCT, 0
731ztitle = "Mean chla"
732PLOT, [ 0.01, 40. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mg chla / m3', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2, TITLE = ztitle
733; add obs
734OPLOT, chla_obs_mean, - depth_obs_mean, PSYM = 8, THICK = 3, SYMSIZE = 0.8
735FOR i = 0, N_dpt - 1 DO BEGIN
736   OPLOT, [ chla_obs_mean(i) - chla_obs_std(i), chla_obs_mean(i) + chla_obs_std(i) ], $
737          [ -depth_obs_mean(i), -depth_obs_mean(i) ]
738ENDFOR
739; plot seawater value
740LOADCT, 1
741OPLOT, [ 0.155 ], [ -1 ], PSYM = 8, THICK = 3, SYMSIZE = 1.6, COLOR = 150
742
743; add model
744LOADCT, 9
745OPLOT, chla_mod_mean, -zz, THICK = 3, COLOR = 150
746LOADCT, 0
747OPLOT, chla_mod_mean+chla_mod_std, -zz, THICK = 1, COLOR = 150
748OPLOT, chla_mod_mean-chla_mod_std, -zz, THICK = 1, COLOR = 150
749
750;--- NOX ---
751
752; prepare plot
753LOADCT, 0
754ztitle = "Mean NOX"
755PLOT, [ 0.00, 30. ], [ -1.0, 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mmmol N / m3', XTICKFORMAT='(I3)', XMINOR = 2, YMINOR = 2, TITLE = ztitle
756; add obs
757OPLOT, nox_obs_mean, - depth_obs_mean, PSYM = 8, THICK = 3, SYMSIZE = 0.8
758FOR i = 0, N_dpt - 1 DO BEGIN
759   OPLOT, [ nox_obs_mean(i) - nox_obs_std(i), nox_obs_mean(i) + nox_obs_std(i) ], $
760          [ -depth_obs_mean(i), -depth_obs_mean(i) ]
761ENDFOR
762
763; plot seawater value
764LOADCT, 1
765OPLOT, [ 26.90 ], [ -1 ], PSYM = 8, THICK = 3, SYMSIZE = 1.6, COLOR = 150
766
767; add model
768LOADCT, 9
769OPLOT, nox_mod_mean, -zz, THICK = 3, COLOR = 150
770LOADCT, 0
771OPLOT, nox_mod_mean+nox_mod_std, -zz, THICK = 1, COLOR = 150
772OPLOT, nox_mod_mean-nox_mod_std, -zz, THICK = 1, COLOR = 150
773
774;--- PO4 ---
775
776; prepare plot
777LOADCT, 0
778ztitle = "Mean PO4"
779PLOT, [ 0.00, 4.0 ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mmmol P / m3', XTICKFORMAT='(I2)', XMINOR = 2, YMINOR = 2, TITLE = ztitle
780; add obs
781OPLOT, po4_obs_mean, - depth_obs_mean, PSYM = 8, THICK = 3, SYMSIZE = 0.8
782FOR i = 0, N_dpt - 1 DO BEGIN
783   OPLOT, [ po4_obs_mean(i) - po4_obs_std(i), po4_obs_mean(i) + po4_obs_std(i) ], $
784          [ -depth_obs_mean(i), -depth_obs_mean(i) ]
785ENDFOR
786
787; plot seawater value
788LOADCT, 1
789OPLOT, [ 2.35 ], [ -1 ], PSYM = 8, THICK = 3, SYMSIZE = 1.6, COLOR = 150
790
791; add model
792LOADCT, 9
793OPLOT, po4_mod_mean, -zz, THICK = 3, COLOR = 150
794LOADCT, 0
795OPLOT, po4_mod_mean+po4_mod_std, -zz, THICK = 1, COLOR = 150
796OPLOT, po4_mod_mean-po4_mod_std, -zz, THICK = 1, COLOR = 150
797
798;--- DSi ---
799
800; prepare plot
801LOADCT, 0
802ztitle = "Mean DSi"
803PLOT, [ 0.00, 60. ], [ -1., 0. ], CHARSIZE = cs, /NODATA, SUBTITLE = 'mmmol N / m3', XTICKFORMAT='(I2)', XMINOR = 2, YMINOR = 2, TITLE = ztitle
804; add obs
805OPLOT, dsi_obs_mean, - depth_obs_mean, PSYM = 8, THICK = 3, SYMSIZE = 0.8
806FOR i = 0, N_dpt - 1 DO BEGIN
807   OPLOT, [ dsi_obs_mean(i) - dsi_obs_std(i), dsi_obs_mean(i) + dsi_obs_std(i) ], $
808          [ -depth_obs_mean(i), -depth_obs_mean(i) ]
809ENDFOR
810
811; plot seawater value
812LOADCT, 1
813OPLOT, [ 51.25 ], [ -1 ], PSYM = 8, THICK = 3, SYMSIZE = 1.6, COLOR = 150
814
815; add model
816LOADCT, 9
817OPLOT, dsi_mod_mean, -zz, THICK = 3, COLOR = 150
818LOADCT, 0
819OPLOT, dsi_mod_mean+dsi_mod_std, -zz, THICK = 1, COLOR = 150
820OPLOT, dsi_mod_mean-dsi_mod_std, -zz, THICK = 1, COLOR = 150
821
822;
823;==============================================================================
824; End of the script
825;==============================================================================
826;
827IF ( device EQ 'PS' ) THEN BEGIN
828   DEVICE, /CLOSE
829   SET_PLOT, "X"
830   !P.MULTI=[0,2,2]
831ENDIF
832
833END
834
Note: See TracBrowser for help on using the repository browser.