| 1 | pro wind_esa |
|---|
| 2 | |
|---|
| 3 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; |
|---|
| 4 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; |
|---|
| 5 | filename='ok_mesoscale_profile' |
|---|
| 6 | retrieve='true' |
|---|
| 7 | ;retrieve='false' |
|---|
| 8 | file='/d6/vblmd/MERIDIANI_EXOMARS/saves_simu_LS_68.5/jour20/wrfout_d01_2024-03-20_06:00:00' |
|---|
| 9 | start_time=6. |
|---|
| 10 | interv=1. |
|---|
| 11 | landing_site=[353.87-360.,-1.88] |
|---|
| 12 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; |
|---|
| 13 | ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; |
|---|
| 14 | |
|---|
| 15 | if (retrieve eq 'true') then begin |
|---|
| 16 | |
|---|
| 17 | ; |
|---|
| 18 | ; load WRF coordinates |
|---|
| 19 | ; |
|---|
| 20 | getcdf, file=file, charvar='XLAT', invar=lat |
|---|
| 21 | getcdf, file=file, charvar='XLONG', invar=lon |
|---|
| 22 | getcdf, file=file, charvar='HGT', invar=elevation |
|---|
| 23 | ; |
|---|
| 24 | ; find simulated lander grid point |
|---|
| 25 | ; |
|---|
| 26 | findxy, $ |
|---|
| 27 | x=lon, $ |
|---|
| 28 | y=lat, $ |
|---|
| 29 | point=landing_site, $ |
|---|
| 30 | ind=indp, $ |
|---|
| 31 | flag=2 |
|---|
| 32 | ; |
|---|
| 33 | ; get field |
|---|
| 34 | ; |
|---|
| 35 | getcdf, file=file, charvar='U', invar=u |
|---|
| 36 | getcdf, file=file, charvar='V', invar=v |
|---|
| 37 | ;;; BEWARE: arrays not the same size... because staggering... u^2+v^2 is wrong and IDL does not complain... |
|---|
| 38 | field=sqrt( u(indp[0],indp[1],*,*)^2 + v(indp[0],indp[1],*,*)^2 ) & u=0. & v=0. |
|---|
| 39 | what_I_plot=reform(field) |
|---|
| 40 | ; |
|---|
| 41 | ; local time |
|---|
| 42 | ; |
|---|
| 43 | localtime = float(start_time) + float(interv)*findgen(n_elements(field(0,0,0,*))) ;+ lon(indp[0],indp[1])/15. |
|---|
| 44 | ; |
|---|
| 45 | ; get model height |
|---|
| 46 | ; |
|---|
| 47 | getcdf, file=file, charvar='PHTOT', invar=ph |
|---|
| 48 | height=reform(ph(indp[0],indp[1],*,*)) |
|---|
| 49 | ;height=total(height,2)/n_elements(height(0,*)) |
|---|
| 50 | height=height/1000./3.72 |
|---|
| 51 | ;heightp=height(0:n_elements(height(*))-2)-elevation(indp[0],indp[1])/1000. |
|---|
| 52 | heightp = height - elevation(indp[0],indp[1]) / 1000. |
|---|
| 53 | ; |
|---|
| 54 | ; save |
|---|
| 55 | ; |
|---|
| 56 | save, what_I_plot, localtime, heightp, filename='yeahyeah' |
|---|
| 57 | |
|---|
| 58 | endif else begin |
|---|
| 59 | |
|---|
| 60 | restore, filename='yeahyeah' |
|---|
| 61 | |
|---|
| 62 | endelse |
|---|
| 63 | |
|---|
| 64 | lctime=[7.,8.,9.,10.,11.] |
|---|
| 65 | for ll=0,n_elements(lctime)-1 do begin |
|---|
| 66 | w = where(localtime eq lctime(ll)) |
|---|
| 67 | print, localtime[w(0)] |
|---|
| 68 | tab1 = reform(heightp(*,w(0))) |
|---|
| 69 | tab2 = reform(what_I_plot(*,w(0))) |
|---|
| 70 | tab1 = tab1[0:n_elements(tab1)-14] |
|---|
| 71 | tab2 = tab2[0:n_elements(tab2)-14] |
|---|
| 72 | make_ascii, filename+'_LT'+string(lctime(ll),'(I0)'), tab1, tab2 |
|---|
| 73 | plot, tab1, tab2 |
|---|
| 74 | endfor |
|---|
| 75 | |
|---|
| 76 | |
|---|
| 77 | end |
|---|
