source: TOOLS/CMIP6_FORCING/AER_STRAT/volc.sh @ 4030

Last change on this file since 4030 was 4030, checked in by oboucher, 6 years ago

Remove script that is not used
volc.sh modified for MR

  • Property svn:executable set to *
File size: 19.6 KB
Line 
1dirpwd=`pwd`
2
3vv='v3'
4
5#--choose the resolution
6#--only LR and zoom works fully for now !
7#lmdz='VVLR'
8#lmdz='VLR_L79'
9#lmdz='VLR'
10#lmdz='LR'
11lmdz='MR'
12#--added for Frederique Cheruy
13#lmdz='zoom_128x89'
14
15dirout='/data/'${USER}'/CMIP6/VOLC/'${lmdz}'_'${vv}'/'
16echo $dirout
17if [ ! -d ${dirout} ] ; then mkdir -p ${dirout} ; fi
18
19cat > process_volc.pro << EOF
20pro regrid
21;--version sans longitude
22;--latitude - level - wavelength - time
23;--sep 2017: add lat and month z variations of model layers
24;--Olivier Boucher
25;
26lmdz='${lmdz}'
27;
28if (lmdz eq 'VVLR') then begin
29dimlonlmdz=48
30dimlatlmdz=37
31dimz=19
32output='/data/${USER}/CMIP6/VOLC/VVLR_${vv}/'
33latfirst=90.
34latinc=-180./float(dimlatlmdz-1)
35latitudelmdz=latfirst+latinc*indgen(dimlatlmdz)
36endif
37
38if (lmdz eq 'VLR') then begin
39dimlonlmdz=96
40dimlatlmdz=96
41dimz=39
42output='/data/${USER}/CMIP6/VOLC/VLR_${vv}/'
43latfirst=90.
44latinc=-180./float(dimlatlmdz-1)
45latitudelmdz=latfirst+latinc*indgen(dimlatlmdz)
46endif
47
48if (lmdz eq 'VLR_L79') then begin
49dimlonlmdz=96
50dimlatlmdz=96
51dimz=79
52output='/data/${USER}/CMIP6/VOLC/VLR_L79_${vv}/'
53latfirst=90.
54latinc=-180./float(dimlatlmdz-1)
55latitudelmdz=latfirst+latinc*indgen(dimlatlmdz)
56endif
57
58if (lmdz eq 'LR') then begin
59dimlonlmdz=144
60dimlatlmdz=143
61dimz=79
62output='/data/${USER}/CMIP6/VOLC/LR_${vv}/'
63latfirst=90.
64latinc=-180./float(dimlatlmdz-1)
65latitudelmdz=latfirst+latinc*indgen(dimlatlmdz)
66endif
67
68if (lmdz eq 'MR') then begin
69dimlonlmdz=256
70dimlatlmdz=257
71dimz=79
72output='/data/${USER}/CMIP6/VOLC/MR_${vv}/'
73latfirst=90.
74latinc=-180./float(dimlatlmdz-1)
75latitudelmdz=latfirst+latinc*indgen(dimlatlmdz)
76endif
77
78if (lmdz eq 'zoom_128x89') then begin
79dimlonlmdz=128
80dimlatlmdz=89
81dimz=79
82output='/data/${USER}/CMIP6/VOLC/${lmdz}_${vv}/'
83latitudelmdz=[90., 86.85495, 83.8430634, 80.9792252, 78.2779922, 75.7526321, 73.4137878, 71.2680206, 69.3164978, $
84              67.5542145,   65.9702301, 64.5487747, 63.2709846, 62.1167107, 61.0660896, 60.1005592, 59.2034454, $
85              58.360157,    57.5581551, 56.7868576, 56.0374718, 55.302887,  54.5775719, 53.8574829, 53.1399384, $
86              52.4233971,   51.7071419, 50.9909439, 50.2747459, 49.5585518, 48.8423538, 48.1261597, 47.4099655, $
87              46.6937675,   45.9775734, 45.2613792, 44.5451813, 43.8289604, 43.1125679, 42.395504,  41.6765518, $
88              40.9534645,   40.2226982, 39.4792557, 38.7165909, 37.9265213, 37.099102,  36.2225075, 35.2828789, $
89              34.2642632,   33.1487236, 31.9167614, 30.5481758, 29.0234509, 27.32551,   25.4415455, 23.3643932, $
90              21.0930786,   18.6323872, 15.9917059, 13.1835403, 10.2221146, 7.12228632, 3.89880991, 0.565921485,$
91              -2.86286402, -6.37483072, -9.95808029, -13.6015415, -17.2949867, -21.029068, -24.7953758, -28.5864811, $
92             -32.3959961, -36.2186012, -40.0500298, -43.8870316, -47.7272758, -51.5692215, -55.4119453, -59.2549667, $
93             -63.0980797, -66.9412079, -70.7843399, -74.6274719, -78.4706039, -82.313736, -86.156868, -90.]
94endif
95
96NSW=6
97NLW=16
98;
99;--script only works on ciclad
100dir='/prodigfs/project/input4MIPs/VOLC/${vv}/'
101filename=dir+'CMIP_1850_2014_extinction_550nm_${vv}.nc'
102NETCDFREAD,filename,'altitude',altitude,dimaltitude
103NETCDFREAD,filename,'ext550',ext550,dimext
104;
105filename=dir+'CMIP_IPSL-CM6_radiation_${vv}.nc'
106;
107NETCDFREAD,filename,'altitude',altitude,dimaltitude
108NETCDFREAD,filename,'latitude',latitude,dimlatitude
109NETCDFREAD,filename,'wl1_sun',wl1_sun,solar_bands
110NETCDFREAD,filename,'wl2_sun',wl2_sun,solar_bands
111NETCDFREAD,filename,'wl1_earth',wl1_earth,terrestrial_bands
112NETCDFREAD,filename,'wl2_earth',wl2_earth,terrestrial_bands
113NETCDFREAD,filename,'month',month,dimmonth
114NETCDFREAD,filename,'ext_sun',ext_sun,dimext
115NETCDFREAD,filename,'omega_sun',omega_sun,dimext
116NETCDFREAD,filename,'g_sun',g_sun,dimext
117NETCDFREAD,filename,'ext_earth',ext_earth,dimext
118NETCDFREAD,filename,'omega_earth',omega_earth,dimext
119NETCDFREAD,filename,'g_earth',g_earth,dimext
120;
121print ,'wl1_sun=', wl1_sun
122print ,'wl2_sun=', wl2_sun
123print ,'wl1_ear=', wl1_earth
124print ,'wl2_ear=', wl2_earth
125;
126print ,'min max ext_sun=', min(ext_sun),   max(ext_sun)
127print ,'min max ome_sun=', min(omega_sun), max(omega_sun)
128;for nl=0, NSW-1 do begin
129;print ,'min max ome_sun nl=', nl, min(omega_sun(*,*,*,nl)), max(omega_sun(*,*,*,nl))
130;endfor
131print ,'min max ggg_sun=', min(g_sun),     max(g_sun)
132print ,'min max ext_ear=', min(ext_earth), max(ext_earth)
133;
134dimlat=dimlatitude(0)
135dimalt=dimaltitude(0)
136dimtime=dimmonth(0)
137;
138;--vertical resolution in Beiping's code
139dz=0.5
140;
141month_in_year=12
142;;month_in_year=1
143;
144;--compute optical depth of input data
145;--dimension ext_sun(l,k,j,nl)
146;ext_sun0=ext_sun(0,*,*,*)
147;print , 'size ext_sun0=', size(ext_sun0)
148;tau_sun=TOTAL(ext_sun0,2)*dz
149;print ,'min max tau_sun l=0 =', min(tau_sun),   max(tau_sun)
150;
151zz=fltarr(dimz)
152;
153;---approximate altitudes of LMDZ - L19 and L39
154;---NOT RECOMMENDED
155;
156if (dimz eq 19) then begin
157zz=[ 0.07125416, 0.2402302, 0.4865242, 0.8627044, 1.426484, 2.244528,    $
158    3.394085, 4.952366, 6.950302, 9.26943, 11.59338, 13.71214, 15.83651, $
159    18.27338, 21.16397, 24.63156, 28.9039, 34.57184, 44.51565] ;
160endif
161;
162if (dimz eq 39) then begin
163zz=[0.0338988, 0.1106602, 0.2139233, 0.3609663, 0.5685416, 0.8534464,     $
164    1.233232, 1.726868, 2.355076, 3.139813, 4.101937, 5.25541, 6.596076,  $
165    8.085103, 9.636482, 11.13441, 12.50023, 13.75765, 15.00424, 16.32207, $
166    17.74403, 19.27899, 20.93219, 22.70941, 24.61736, 26.66389, 28.85831, $
167    31.21178, 33.73778, 36.45278, 39.37709, 42.53607, 45.96209, 49.69786, $
168    53.80379, 58.37727, 63.61497, 70.07092, 80.52708]
169endif
170;
171;---exact altitudes of LMDZ -- L79
172;--only exists for LR, MR and zoom_128x89 at the moment
173if (dimz eq 79) then begin
174  filename='./zalt_zonmean_${lmdz}_l79_rev.nc'
175  NETCDFREAD,filename,'GEOP',zz,dimzz
176  NETCDFREAD,filename,'LAT',zzlat,dimzzlat
177  NETCDFREAD,filename,'TIME_COUNTER',zztime,dimzztime
178  dimzzilat=dimzzlat[0]
179  dimzzitime=dimzztime[0]
180  ;--becareful zz comes with four dimensions
181  ;--lon lat k time
182  print, 'GEOP size=', size(zz)
183  ;--becareful zzlat comes with South Pole first
184  print, 'LAT from zalt field=', zzlat
185  if (dimzzilat ne dimlatlmdz) then begin
186    print , 'PB dimension latitude'
187  endif
188endif
189;
190;
191;--reconstructing the vertical coordinate at interfaces (in unit km)
192;--and reverse lat axis for zzi  i==>ii
193;--and reverse zz axis for zzi   dimz-1-k==>k
194;--and forget about lon axis : index 0 in zz
195zzi=fltarr(dimzzilat,dimzzitime,dimz+1)
196for i=0, dimzzilat-1 do begin
197 ii=dimzzilat-1-i
198 for t=0, dimzzitime-1 do begin
199  zzi(ii,t,0)=0.0
200  for k=1, dimz-1 do begin
201   zzi(ii,t,k)=(zz(0,i,dimz-1-(k-1),t)+zz(0,i,dimz-1-k,t))/2.0
202  endfor
203  zzi(ii,t,dimz)=100.00
204 endfor
205endfor
206print, 'zzi=',zzi
207;
208lev=indgen(dimz)+1
209;
210dimzori=dimalt
211zziori=fltarr(dimzori+1)
212zziori(0)=altitude(0)-0.25
213for k=1, dimzori do begin
214 zziori(k)=altitude(k-1)+0.25
215endfor
216;
217;--550 nm properties
218tau_550_lmdz=fltarr(dimlatlmdz,dimz,month_in_year)
219;
220;--SW properties
221tau_sun_lmdz=fltarr(dimlatlmdz,dimz,NSW,month_in_year)
222ome_sun_lmdz=fltarr(dimlatlmdz,dimz,NSW,month_in_year)
223ggg_sun_lmdz=fltarr(dimlatlmdz,dimz,NSW,month_in_year)
224;
225tau_sun_lmdz_ave=fltarr(dimlatlmdz,dimz,NSW,month_in_year)
226ome_sun_lmdz_ave=fltarr(dimlatlmdz,dimz,NSW,month_in_year)
227ggg_sun_lmdz_ave=fltarr(dimlatlmdz,dimz,NSW,month_in_year)
228;
229tau_sun_lmdz_tr=fltarr(dimlatlmdz,dimz,NSW,month_in_year)
230ome_sun_lmdz_tr=fltarr(dimlatlmdz,dimz,NSW,month_in_year)
231ggg_sun_lmdz_tr=fltarr(dimlatlmdz,dimz,NSW,month_in_year)
232;
233;--LW properties tau_abs
234tau_ear_lmdz=fltarr(dimlatlmdz,dimz,NLW,month_in_year)
235;
236tau_ear_lmdz_ave=fltarr(dimlatlmdz,dimz,NLW,month_in_year)
237;
238tau_ear_lmdz_tr=fltarr(dimlatlmdz,dimz,NLW,month_in_year)
239;
240tau_sun_lmdz_ave(*,*,*,0)=1.e-15
241ome_sun_lmdz_ave(*,*,*,0)=1.e-15
242ggg_sun_lmdz_ave(*,*,*,0)=0.0
243tau_ear_lmdz_ave(*,*,*,0)=1.e-15
244;
245for year=0, dimtime/month_in_year-1 do begin
246;
247tau_550_lmdz(*,*,*)=1.e-15
248tau_sun_lmdz(*,*,*,*)=1.e-15
249ome_sun_lmdz(*,*,*,*)=1.e-15
250ggg_sun_lmdz(*,*,*,*)=0.0
251tau_ear_lmdz(*,*,*,*)=1.e-15
252;
253chyr=strcompress(1850+year,/rem)
254;
255for mth=0,month_in_year-1 do begin
256;
257;--timestep
258l=mth+month_in_year*year
259print,'year mth l=',chyr, mth, l
260;
261;regridding
262for j=0, dimlatlmdz-1  do begin
263;
264;--finding latitude in beiping luo data
265jj=0
266for jluo=0,dimlat-2 do begin
267  if (latitudelmdz(j) gt (latitude(jluo)+latitude(jluo+1))/2. ) then begin
268    jj=jluo+1
269  endif
270endfor
271;
272for k=0, dimz-1 do begin
273for kori=0, dimzori-1 do begin
274;
275;fraction de la maille kori qui se trouve dans la maille k
276frac= max([0.0,min([zzi(j,mth,k+1),zziori(kori+1)])-max([zzi(j,mth,k),zziori(kori)])])/(zziori(kori+1)-zziori(kori))
277;
278tau_550_lmdz(j,k,mth)=tau_550_lmdz(j,k,mth)+ext550(l,kori,jj)*dz*frac
279;
280for nl=0, NSW-1 do begin
281  tau_sun_lmdz(j,k,nl,mth)=tau_sun_lmdz(j,k,nl,mth)+ext_sun(l,kori,jj,nl)*dz*frac
282  ome_sun_lmdz(j,k,nl,mth)=ome_sun_lmdz(j,k,nl,mth)+omega_sun(l,kori,jj,nl)*ext_sun(l,kori,jj,nl)*dz*frac
283  ggg_sun_lmdz(j,k,nl,mth)=ggg_sun_lmdz(j,k,nl,mth)+g_sun(l,kori,jj,nl)*omega_sun(l,kori,jj,nl)*ext_sun(l,kori,jj,nl)*dz*frac
284;
285  tau_sun_lmdz_ave(j,k,nl,0)=tau_sun_lmdz_ave(j,k,nl,0)+ext_sun(l,kori,jj,nl)*dz*frac
286  ome_sun_lmdz_ave(j,k,nl,0)=ome_sun_lmdz_ave(j,k,nl,0)+omega_sun(l,kori,jj,nl)*ext_sun(l,kori,jj,nl)*dz*frac
287  ggg_sun_lmdz_ave(j,k,nl,0)=ggg_sun_lmdz_ave(j,k,nl,0)+g_sun(l,kori,jj,nl)*omega_sun(l,kori,jj,nl)*ext_sun(l,kori,jj,nl)*dz*frac
288endfor
289;
290for nl=0, NLW-1 do begin
291  tau_ear_lmdz(j,k,nl,mth)=tau_ear_lmdz(j,k,nl,mth)+ext_earth(l,kori,jj,nl)*(1.-omega_earth(l,kori,jj,nl))*dz*frac
292  tau_ear_lmdz_ave(j,k,nl,0)=tau_ear_lmdz_ave(j,k,nl,0)+ext_earth(l,kori,jj,nl)*(1.-omega_earth(l,kori,jj,nl))*dz*frac
293endfor
294;
295endfor
296;--end lat loop
297;
298endfor
299endfor
300;--end k loops
301;
302endfor
303;--end month loop
304;
305;renormalizing intensive SW properties
306ggg_sun_lmdz(*,*,*,*)=ggg_sun_lmdz(*,*,*,*)/ome_sun_lmdz(*,*,*,*)
307ome_sun_lmdz(*,*,*,*)=ome_sun_lmdz(*,*,*,*)/tau_sun_lmdz(*,*,*,*)
308;
309;saving netcdf file
310;
311print ,'min max ext_sun lmdz=', min(tau_sun_lmdz), max(tau_sun_lmdz)
312print ,'min max ome_sun lmdz=', min(ome_sun_lmdz), max(ome_sun_lmdz)
313;for nl=0, NSW-1 do begin
314;print ,'min max ome_sun lmdz nl =', nl, min(ome_sun_lmdz(*,*,nl,*)), max(ome_sun_lmdz(*,*,nl,*))
315;endfor
316print ,'min max ggg_sun lmdz=', min(ggg_sun_lmdz), max(ggg_sun_lmdz)
317print ,'min max ext_ear lmdz=', min(tau_ear_lmdz), max(tau_ear_lmdz)
318;
319print ,'min max ext_sun lmdz l=1=', min(tau_sun_lmdz(*,*,*,0)), max(tau_sun_lmdz(*,*,*,0))
320print ,'min max ome_sun lmdz l=1=', min(ome_sun_lmdz(*,*,*,0)), max(ome_sun_lmdz(*,*,*,0))
321print ,'min max ggg_sun lmdz l=1=', min(ggg_sun_lmdz(*,*,*,0)), max(ggg_sun_lmdz(*,*,*,0))
322print ,'min max ext_ear lmdz l=1=', min(tau_ear_lmdz(*,*,*,0)), max(tau_ear_lmdz(*,*,*,0))
323;
324;--compute optical depth of output data
325;tau_sun=TOTAL(tau_sun_lmdz,2)
326;print ,'min max tau_sun_lmdz vert=', min(tau_sun),   max(tau_sun)
327;
328opticstruct={lat:fltarr(dimlatlmdz),lev:fltarr(dimz),             $
329             wav:fltarr(NSW),time:fltarr(month_in_year),          $
330             tau_sun:fltarr(dimlatlmdz,dimz,NSW,month_in_year),   $
331             ome_sun:fltarr(dimlatlmdz,dimz,NSW,month_in_year),   $
332             ggg_sun:fltarr(dimlatlmdz,dimz,NSW,month_in_year) }
333;
334opticstruct.lat=latitudelmdz
335opticstruct.lev=lev
336opticstruct.wav=(wl1_sun(0:NSW-1)+wl2_sun(0:NSW-1))/2.
337opticstruct.time=float(indgen(month_in_year)+1)
338opticstruct.tau_sun=tau_sun_lmdz
339opticstruct.ome_sun=ome_sun_lmdz
340opticstruct.ggg_sun=ggg_sun_lmdz
341;
342attributes = {units:strarr(7),long_name:strarr(7)}
343attributes.units =     ['degrees_north','level','meters','month','-','-','-']
344attributes.long_name = ['latitude','level','wavelength','time','tau_sun','ome_sun','g_sun']
345;
346dimensions = {isdim:intarr(7), links:intarr(4,7)}
347       dimensions.isdim =  [1,1,1,1,0,0,0]  ; (1=dimension, 0=variable)
348       dimensions.links = [ [-1,-1,-1,-1],[-1,-1,-1,-1],   $
349                            [-1,-1,-1,-1],[-1,-1,-1,-1],   $
350                            [0,1,2,3],[0,1,2,3],[0,1,2,3]  ]
351;
352netcdfwrite,output+'tauswstrat.2D.'+chyr+'.nc',opticstruct,clobber=1, $
353            attributes=attributes, dimensions=dimensions
354;
355opticstruct={lat:fltarr(dimlatlmdz),lev:fltarr(dimz),             $
356             wav:fltarr(NLW),time:fltarr(month_in_year),          $
357             tau_ear:fltarr(dimlatlmdz,dimz,NLW,month_in_year) }
358;
359opticstruct.lat=latitudelmdz
360opticstruct.lev=lev
361opticstruct.wav=(wl1_earth+wl2_earth)/2.
362opticstruct.time=float(indgen(month_in_year)+1)
363opticstruct.tau_ear=tau_ear_lmdz
364;
365attributes = {units:strarr(5),long_name:strarr(5)}
366attributes.units = ['degrees_north','level','cm-1','month','-']
367attributes.long_name = ['latitude','level','wavenumber','time','tau_ear']
368;
369dimensions = {isdim:intarr(5), links:intarr(4,5)}
370       dimensions.isdim =  [1,1,1,1,0]  ; (1=dimension, 0=variable)
371       dimensions.links = [ [-1,-1,-1,-1],[-1,-1,-1,-1],   $
372                            [-1,-1,-1,-1],[-1,-1,-1,-1],   $
373                            [0,1,2,3]  ]
374;
375netcdfwrite,output+'taulwstrat.2D.'+chyr+'.nc',opticstruct,clobber=1, $
376            attributes=attributes, dimensions=dimensions
377;
378opticstruct={lat:fltarr(dimlatlmdz),lev:fltarr(dimz),             $
379             time:fltarr(month_in_year),                          $
380             tau550:fltarr(dimlatlmdz,dimz,month_in_year) }
381;
382opticstruct.lat=latitudelmdz
383opticstruct.lev=lev
384opticstruct.time=float(indgen(month_in_year)+1)
385opticstruct.tau550=tau_550_lmdz
386;
387attributes = {units:strarr(4),long_name:strarr(4)}
388attributes.units = ['degrees_north','level','month','-']
389attributes.long_name = ['latitude','level','time','tau550']
390;
391dimensions = {isdim:intarr(4), links:intarr(3,4)}
392       dimensions.isdim =  [1,1,1,0]  ; (1=dimension, 0=variable)
393       dimensions.links = [ [-1,-1,-1],[-1,-1,-1],   $
394                            [-1,-1,-1],[0,1,2]  ]
395;
396netcdfwrite,output+'tau550strat.2D.'+chyr+'.nc',opticstruct,clobber=1, $
397            attributes=attributes, dimensions=dimensions
398;
399endfor
400;--end loop on years
401;
402;now deal with average conditions
403ggg_sun_lmdz_ave(*,*,*,0)=ggg_sun_lmdz_ave(*,*,*,0)/ome_sun_lmdz_ave(*,*,*,0)
404ome_sun_lmdz_ave(*,*,*,0)=ome_sun_lmdz_ave(*,*,*,0)/tau_sun_lmdz_ave(*,*,*,0)
405tau_sun_lmdz_ave(*,*,*,0)=tau_sun_lmdz_ave(*,*,*,0)/float(dimtime)
406tau_ear_lmdz_ave(*,*,*,0)=tau_ear_lmdz_ave(*,*,*,0)/float(dimtime)
407;
408for mth=1, month_in_year-1 do begin
409  ggg_sun_lmdz_ave(*,*,*,mth)=ggg_sun_lmdz_ave(*,*,*,0)
410  ome_sun_lmdz_ave(*,*,*,mth)=ome_sun_lmdz_ave(*,*,*,0)
411  tau_sun_lmdz_ave(*,*,*,mth)=tau_sun_lmdz_ave(*,*,*,0)
412  tau_ear_lmdz_ave(*,*,*,mth)=tau_ear_lmdz_ave(*,*,*,0)
413endfor
414;
415opticstruct={lat:fltarr(dimlatlmdz),lev:fltarr(dimz),             $
416             wav:fltarr(NSW),time:fltarr(month_in_year),          $
417             tau_sun:fltarr(dimlatlmdz,dimz,NSW,month_in_year),   $
418             ome_sun:fltarr(dimlatlmdz,dimz,NSW,month_in_year),   $
419             ggg_sun:fltarr(dimlatlmdz,dimz,NSW,month_in_year) }
420;
421opticstruct.lat=latitudelmdz
422opticstruct.lev=lev
423opticstruct.wav=(wl1_sun(0:NSW-1)+wl2_sun(0:NSW-1))/2.
424opticstruct.time=float(indgen(month_in_year))
425opticstruct.tau_sun=tau_sun_lmdz_ave
426opticstruct.ome_sun=ome_sun_lmdz_ave
427opticstruct.ggg_sun=ggg_sun_lmdz_ave
428;
429attributes = {units:strarr(7),long_name:strarr(7)}
430attributes.units =     ['degrees_north','level','meters','month','-','-','-']
431attributes.long_name = ['latitude','level','wavelength','time','tau_sun','ome_sun','g_sun']
432;
433dimensions = {isdim:intarr(7), links:intarr(4,7)}
434       dimensions.isdim =  [1,1,1,1,0,0,0]  ; (1=dimension, 0=variable)
435       dimensions.links = [ [-1,-1,-1,-1],[-1,-1,-1,-1],   $
436                            [-1,-1,-1,-1],[-1,-1,-1,-1],   $
437                            [0,1,2,3],[0,1,2,3],[0,1,2,3]  ]
438;
439netcdfwrite,output+'tauswstrat.2D.ave.nc',opticstruct,clobber=1, $
440            attributes=attributes, dimensions=dimensions
441;
442opticstruct={lat:fltarr(dimlatlmdz),lev:fltarr(dimz),             $
443             wav:fltarr(NLW),time:fltarr(month_in_year),          $
444             tau_ear:fltarr(dimlatlmdz,dimz,NLW,month_in_year) }
445;
446opticstruct.lat=latitudelmdz
447opticstruct.lev=lev
448opticstruct.wav=(wl1_earth+wl2_earth)/2.
449opticstruct.time=float(indgen(month_in_year))
450opticstruct.tau_ear=tau_ear_lmdz_ave
451;
452attributes = {units:strarr(5),long_name:strarr(5)}
453attributes.units = ['degrees_north','level','cm-1','month','-']
454attributes.long_name = ['latitude','level','wavenumber','time','tau_ear']
455;
456dimensions = {isdim:intarr(5), links:intarr(4,5)}
457       dimensions.isdim =  [1,1,1,1,0]  ; (1=dimension, 0=variable)
458       dimensions.links = [ [-1,-1,-1,-1],[-1,-1,-1,-1],   $
459                            [-1,-1,-1,-1],[-1,-1,-1,-1],   $
460                            [0,1,2,3]  ]
461;
462netcdfwrite,output+'taulwstrat.2D.ave.nc',opticstruct,clobber=1, $
463            attributes=attributes, dimensions=dimensions
464;
465;--finally prepare the transition period 2015-2023
466;--we mix the last year of the period with the climatology
467;
468for year=2015,2023 do begin
469;
470chyr=strcompress(year,/rem)
471;--compute weights
472w1=float(2024-year)/10.
473w2=1.0-w1
474;--SW properties
475tau_sun_lmdz_tr=w1*tau_sun_lmdz+w2*tau_sun_lmdz_ave
476ome_sun_lmdz_tr=w1*tau_sun_lmdz*ome_sun_lmdz+w2*tau_sun_lmdz_ave*ome_sun_lmdz_ave
477ggg_sun_lmdz_tr=w1*tau_sun_lmdz*ome_sun_lmdz*ggg_sun_lmdz+ $
478                w2*tau_sun_lmdz_ave*ome_sun_lmdz_ave*ggg_sun_lmdz_ave
479ggg_sun_lmdz_tr=ggg_sun_lmdz_tr/ome_sun_lmdz_tr
480ome_sun_lmdz_tr=ome_sun_lmdz_tr/tau_sun_lmdz_tr
481;
482;--LW properties
483tau_ear_lmdz_tr=w1*tau_ear_lmdz+w2*tau_ear_lmdz_ave
484;
485;--prepare SW output
486opticstruct={lat:fltarr(dimlatlmdz),lev:fltarr(dimz),             $
487             wav:fltarr(NSW),time:fltarr(month_in_year),          $
488             tau_sun:fltarr(dimlatlmdz,dimz,NSW,month_in_year),   $
489             ome_sun:fltarr(dimlatlmdz,dimz,NSW,month_in_year),   $
490             ggg_sun:fltarr(dimlatlmdz,dimz,NSW,month_in_year) }
491;
492opticstruct.lat=latitudelmdz
493opticstruct.lev=lev
494opticstruct.wav=(wl1_sun(0:NSW-1)+wl2_sun(0:NSW-1))/2.
495opticstruct.time=float(indgen(month_in_year)+1)
496opticstruct.tau_sun=tau_sun_lmdz_tr
497opticstruct.ome_sun=ome_sun_lmdz_tr
498opticstruct.ggg_sun=ggg_sun_lmdz_tr
499;
500attributes = {units:strarr(7),long_name:strarr(7)}
501attributes.units =     ['degrees_north','level','meters','month','-','-','-']
502attributes.long_name = ['latitude','level','wavelength','time','tau_sun','ome_sun','g_sun']
503;
504dimensions = {isdim:intarr(7), links:intarr(4,7)}
505       dimensions.isdim =  [1,1,1,1,0,0,0]  ; (1=dimension, 0=variable)
506       dimensions.links = [ [-1,-1,-1,-1],[-1,-1,-1,-1],   $
507                            [-1,-1,-1,-1],[-1,-1,-1,-1],   $
508                            [0,1,2,3],[0,1,2,3],[0,1,2,3]  ]
509;
510netcdfwrite,output+'tauswstrat.2D.'+chyr+'.nc',opticstruct,clobber=1, $
511            attributes=attributes, dimensions=dimensions
512;
513;--prepare LW output
514opticstruct={lat:fltarr(dimlatlmdz),lev:fltarr(dimz),             $
515             wav:fltarr(NLW),time:fltarr(month_in_year),          $
516             tau_ear:fltarr(dimlatlmdz,dimz,NLW,month_in_year) }
517;
518opticstruct.lat=latitudelmdz
519opticstruct.lev=lev
520opticstruct.wav=(wl1_earth+wl2_earth)/2.
521opticstruct.time=float(indgen(month_in_year)+1)
522opticstruct.tau_ear=tau_ear_lmdz_tr
523;
524attributes = {units:strarr(5),long_name:strarr(5)}
525attributes.units = ['degrees_north','level','cm-1','month','-']
526attributes.long_name = ['latitude','level','wavenumber','time','tau_ear']
527;
528dimensions = {isdim:intarr(5), links:intarr(4,5)}
529       dimensions.isdim =  [1,1,1,1,0]  ; (1=dimension, 0=variable)
530       dimensions.links = [ [-1,-1,-1,-1],[-1,-1,-1,-1],   $
531                            [-1,-1,-1,-1],[-1,-1,-1,-1],   $
532                            [0,1,2,3]  ]
533;
534netcdfwrite,output+'taulwstrat.2D.'+chyr+'.nc',opticstruct,clobber=1, $
535            attributes=attributes, dimensions=dimensions
536;
537endfor
538;--end loop on years
539;
540end
541EOF
542
543cat > volc.job << EOF2
544#PBS -N process_volc
545#PBS -S /bin/bash
546#PBS -q week # there exist: short, day, days3, week...
547#PBS -k eo  # to write the output of stdin
548### Max memory
549#PBS -l vmem=10gb   # virtual memory
550#PBS -l  mem=10gb
551
552cd $dirpwd
553
554idl << eof
555.r netcdf
556.r process_volc
557regrid
558eof
559EOF2
560
561qsub volc.job
Note: See TracBrowser for help on using the repository browser.