source: trunk/forcageforMOED.pro @ 2

PRO forcageforMOED
@init2
@initorca2_bab

@common

ian='01'
iyear='1993'

e_exp='ESS'
rep_fred='/usr/work/sur/fvi/OPA/geomag/'
key_portrait = 0
; stockage des fichiers brut 
  ioDATA='/usr/work/sur/fvi/OPA/ORCA2/'
  file_U=e_exp+'_5d_'+ian+'0101_'+ian+'1231_grid_U.nc' 
  print, file_U
  file_V=e_exp+'_5d_'+ian+'0101_'+ian+'1231_grid_V.nc' 
  print, file_V
  file_T=e_exp+'_5d_'+ian+'0101_'+ian+'1231_grid_T.nc' 
  print, file_T
  file_Sed= rep_fred+'cond_sed_ORCA2.nc'
  file_Br= rep_fred+'Br_ORCA2.nc'


; title
     t_exp= e_exp
        
     t_bt  = 'bar_transp'
ioORLN2 = '/usr/work/sur/fvi/OPA/ORCA2'
;facteur d'echelle vertical  for partial steps
e3v3d=read_ncdf('e3v_ps',0,/timestep,iodir=ioORLN2,/nostruct,/tout,filename='meshmask_bab.nc')
e3u3d=read_ncdf('e3u_ps',0,/timestep,iodir=ioORLN2,/nostruct,/tout,filename='meshmask_bab.nc')
SIGMAsed=read_ncdf('cond_sed',0,/timestep,/nostruct,/tout,filename=file_Sed,/cont_nofill)
BR=read_ncdf('Br',0,/timestep,/nostruct,/tout,filename=file_Br,/cont_nofill)





; vertical integration:
e3t3d=make_array(jpi,jpj,jpk)   
        for k=0, jpk-1 do begin              &$
          for j=0,jpj-1 do begin              &$      
            for i=0,jpi-1 do begin             &$
              e3t3d(i,j,k) = e3t(k)    &$
            endfor                     &$
          endfor                      &$
        endfor
jpt = 73

;vud = make_array(jpi,jpj,jpt)   
;vvd = make_array(jpi, jpj, jpt)
divBustar = make_array(jpi, jpj, jpt)
diver3=replicate(0.,182,149,1,73)
sigma3=replicate(0.,182,149,1,73)
Jpu=replicate(0.,182,149,1,73)
Jpv=replicate(0.,182,149,1,73)


FOR jt = 0, jpt-1 DO BEGIN &$   

; ouverture des fichiers et stockage en memoire partial steps
  vu=read_ncdf('vozocrtx',jt,jt, /timestep, iodir=ioDATA,/nostruct,/TOUT,filename=file_U)  &$   
;stop
  vv=read_ncdf('vomecrty',jt,jt, /timestep,iodir=ioDATA,/nostruct,/TOUT,filename=file_V)  &$   
;stop
; lecture salinite & temperature
  temp= read_ncdf('votemper',jt,jt,/timestep,iodir=ioDATA,/nostruct,/tout,filename=file_T)
;stop
  salin=read_ncdf('vosaline',jt,jt,/timestep,iodir=ioDATA,/nostruct,/tout,filename=file_T)
;stop
  conduct=0.02047780622061 + 0.00273147624197*temp + 0.00035133182334*temp*temp + 0.09139808809909*salin + 0.00241425798890*salin*temp -0.00023998958774*salin*salin
  mask_t=where(conduct GT 1.e+19)
  conduct(mask_t)=0.
; Somme conduct au point T


; Calcul SIGMA

  SIGMAoc=total(conduct*e3t3d*tmask,3)
  SIGMA=SIGMAsed+SIGMAoc



  SIGMA_u=(SIGMA+shift(SIGMA,-1,0))/2.

  SIGMA_v=(SIGMA+shift(SIGMA,0,1))/2.

; Calcul B en points u et v

  BR_u=(BR+shift(BR,-1,0))/2.

  BR_v=(BR+shift(BR,0,1))/2.


; Calcul integrale conduct

  conduct_u=(conduct+shift(conduct,-1,0,0))/2.

  conduct_v=(conduct+shift(conduct,0,1,0))/2.
  
  u_cond_u=total( vu*conduct_u*e3u3d*umask(),3)

  v_cond_v=total( vv*conduct_v*e3v3d*vmask(),3)


;calcul de Jp= u_star x Fz k
  Jpv_star= -1*BR_u*u_cond_u

  Jpu_star= BR_v*v_cond_v

; Divergence du champ

  Diver=divfred(Jpu_star,Jpv_star)

Diver=Diver*1e-6
;stop
lecontinent=where(Diver GT 1.E08)
Diver(lecontinent)=0.



;stop
;bande de recouvrement::

diver3(1:180,0:147,*,jt)=Diver(*,*)
diver3(0,*,*,*)=diver3(180,*,*,*)
diver3(181,*,*,*)=diver3(1,*,*,*)

sigma3(1:180,0:147,*,jt)=SIGMA(*,*)
sigma3(0,*,*,*)=sigma3(180,*,*,*)
sigma3(181,*,*,*)=sigma3(1,*,*,*)

Jpu(1:180,0:147,*,jt)=Jpu_star(*,*)
Jpu(0,*,*,*)=Jpu(180,*,*,*)
Jpu(181,*,*,*)=Jpu(1,*,*,*)

Jpv(1:180,0:147,*,jt)=Jpv_star(*,*)
Jpv(0,*,*,*)=Jpv(180,*,*,*)
Jpv(181,*,*,*)=Jpv(1,*,*,*)




print,  jt


ENDFOR 
 



temps=fltarr(73)
temps(0)=0.
for jt=0,71 do begin &$
temps(jt+1)=temps(jt) +5*86400. &$
endfor
print,temps

   vargrid = 'T'
   iodir = '/usr/work/sur/fvi/OPA/ORCA2/'
; Nom
   idout = NCDF_CREATE(iodir+'ForcageMOED_5d_'+iyear+'_grid_T.nc',/clobber)
   print, 'Creation du fichier Netcdf'
   NCDF_CONTROL, idout, /nofill
; Dimension
   xidout = NCDF_DIMDEF(idout, 'x',jpiglo)
   yidout = NCDF_DIMDEF(idout, 'y',jpjglo)
   didout = NCDF_DIMDEF(idout, 'deptht',1)
   tidout = NCDF_DIMDEF(idout, 'time_counter', /unlimited)

   didout1 = NCDF_DIMDEF(idout, 'deptht1',jpk)


; Attributs globaux
   id0  = NCDF_VARDEF(idout, 'nav_lon'     , [xidout, yidout                ], /FLOAT)
   id1  = NCDF_VARDEF(idout, 'nav_lat'     , [xidout, yidout                ], /FLOAT)
   id2  = NCDF_VARDEF(idout, 'deptht'      , [                didout1        ], /FLOAT)
   id3  = NCDF_VARDEF(idout, 'time_counter', [                        tidout], /FLOAT)
   id4  = NCDF_VARDEF(idout, 'DivJp'  , [xidout, yidout, didout, tidout], /DOUBLE)
   id5  = NCDF_VARDEF(idout, 'Sigma'  , [xidout, yidout, didout, tidout], /DOUBLE)
   id6  = NCDF_VARDEF(idout, 'Jpu'  , [xidout, yidout, didout, tidout], /DOUBLE)
   id7  = NCDF_VARDEF(idout, 'Jpv'  , [xidout, yidout, didout, tidout], /DOUBLE)

; Variable 0
   NCDF_ATTPUT, idout, id0, 'units', 'degrees_east'
   NCDF_ATTPUT, idout, id0, 'long_name', 'Longitude'
   NCDF_ATTPUT, idout, id0, 'nav_model', 'Default grid'
; Variable 1
   NCDF_ATTPUT, idout, id1, 'units', 'degrees_north'
   NCDF_ATTPUT, idout, id1, 'long_name', 'Latitude'
   NCDF_ATTPUT, idout, id1, 'nav_model', 'Default grid'
; Variable 2
   NCDF_ATTPUT, idout, id2, 'units','meters'
   NCDF_ATTPUT, idout, id2, 'long_name','Depth'
   NCDF_ATTPUT, idout, id2, 'nav_model','Default grid'
; Variable3
   NCDF_ATTPUT, idout, id3, 'units', 'seconds since 0001-01-01 00:00:00 '
   NCDF_ATTPUT, idout, id3, 'calendar','noleap'
   NCDF_ATTPUT, idout, id3, 'title', 'Time'
   NCDF_ATTPUT, idout, id3, 'long_name', 'Time axis'
   NCDF_ATTPUT, idout, id3, 'time_origin','0001-JAN-01 00:00:00'
; Variables
   NCDF_ATTPUT, idout, id4, 'long_name', 'Divergence Jp'
; Variables
   NCDF_ATTPUT, idout, id5, 'long_name', 'Total Conductance'
   NCDF_ATTPUT, idout, id6, 'long_name', 'Jpu=u_star x B (x comp)'
   NCDF_ATTPUT, idout, id7, 'long_name', 'Jpv=u_star x B (y comp)'


   NCDF_CONTROL, idout, /ENDEF

; Creation de la longitude
   NCDF_VARPUT, idout, id0, glamt
; Creation de la latitude
   NCDF_VARPUT, idout, id1, gphit
; Creation de la profondeur 
   NCDF_VARPUT, idout, id2, gdept 
; Creation du calendrier

   NCDF_VARPUT, idout, id3, temps

   
; Ecriture des donnees

; ecriture des glam_8 
   NCDF_VARPUT, idout, id4 , diver3
   NCDF_VARPUT, idout, id5 , sigma3
   NCDF_VARPUT, idout, id6 , Jpu
   NCDF_VARPUT, idout, id7 , Jpv



   NCDF_CLOSE, idout


END