Changeset 173 for trunk/tools


Ignore:
Timestamp:
11/18/09 10:35:21 (15 years ago)
Author:
pinsard
Message:

split files with multiple pro or function (except wavelets because external)

Location:
trunk/tools
Files:
4 edited
4 copied

Legend:

Unmodified
Added
Removed

  • trunk/tools/bsf/build_laplacien_f.pro

    r170 r173  
    149149      ENDIF  
    150150 
    151       colonne = ligne+decalage[ordre(i MOD 4)-1] 
     151      colonne = ligne+decalage[ordre[i MOD 4]-1] 
    152152; 
    153153      IF (colonne GE 0) AND (colonne LT N2) THEN  BEGIN  

  • trunk/tools/bsf/diag_bsf.pro

    r170 r173  
    11;+ 
    22; 
    3 ; Calcul de la Fonction de Courant Barotrope  
    4 ; a partir du champ de vitesse (u,v) 
     3; @param u 
     4; 
     5; @param v 
     6; 
     7; @keyword FLUX 
     8; 
     9; @keyword _EXTRA 
    510; 
    611; @uses 
    712; <pro>common</pro> 
    8 ; <propost_it>com_eg</propost_it> 
     13; <propost_it>bsf_lec</propost_it> 
    914; 
    1015; @history 
     16; - fplod 20091118T075330Z  
     17; 
     18;   * externalize diag_bsf_mu and diag_bsf_lec 
     19; 
    1120; creation : printemps 98 G. Roullet 
    1221; 
    1322; @version 
    1423; $Id$ 
    15 ; 
    16 ;- 
    17 PRO diag_bsf_mu, hu, hv, mu 
    18 ; 
    19   compile_opt idl2, strictarrsubs 
    20 ; 
    21 @common 
    22 @com_eg 
    23 ; 
    24 ;; calcul adapte a la grille ORCA 
    25 ;; 
    26    mu = fltarr(14) 
    27 ; 0/ Antartic 
    28 ;   mu[0] = total( (hu*e2u*umaskr[*, *, 0])[149, 5:35]) 
    29    mu[0] = total( (hu*e2u*umaskr[*, *, 0])[151, 14:48]) 
    30 ; 1/ America 
    31 ;   mu[1] = -total( (hv*e1v*vmaskr[*, *, 0])[110:150, 60]) 
    32    mu[1] = -total( (hv*e1v*vmaskr[*, *, 0])[110:150, 60]) 
    33 ; 3/ Australia 
    34 ;   mu[3] = total( (hu*e2u*umaskr[*, *, 0])[17, 40:49]) 
    35    mu[3] = total( (hu*e2u*umaskr[*, *, 0])[20, 50:61]) 
    36 ; 2/ New Zealand 
    37 ;   mu[2] = mu[3]+ total( (hv*e1v*vmaskr[*, *, 0])[30:40, 25]) 
    38    mu[2] = mu[3]+ total( (hv*e1v*vmaskr[*, *, 0])[35:49, 44]) 
    39 ; 4/ Madagascar 
    40 ;   mu[4] = total( (hv*e1v*vmaskr[*, *, 0])[150:159, 40]) 
    41    mu[4] = total( (hv*e1v*vmaskr[*, *, 0])[160:164, 56]) 
    42 ; 5/ New Guinea 
    43 ;   mu[5] = - total( (hu*e2u*umaskr[*, *, 0])[16, 50:55]) 
    44    mu[5] = total( (hv*e2v*vmaskr[*, *, 0])[23:33, 61]) 
    45 ; 6/ Celebes 
    46    mu[6] = total( (hv*e1v*vmaskr[*, *, 0])[13:22, 64]) 
    47 ; 7/ Borneo 
    48 ;   mu[6] = total( (hv*e1v*vmaskr[*, *, 0])[7:12, 60]) 
    49    mu[7] = total( (hv*e1v*vmaskr[*, *, 0])[14:17, 67]) 
    50 ; 8/ Philippines 
    51 ;   mu[7] = total( (hv*e1v*vmaskr[*, *, 0])[10:16, 90]) 
    52    mu[8] = total( (hv*e1v*vmaskr[*, *, 0])[14:22, 86]) 
    53 ; 9/ Tasmania 
    54    mu[9] = mu[3]+ total( (hu*e1u*umaskr[*, *, 0])[34, 41:45]) 
    55 ; 10/ Cuba 
    56 ;   mu[8] = - total( (hv*e1v*vmaskr[*, *, 0])[101:127, 86]) 
    57    mu[10] = - total( (hv*e1v*vmaskr[*, *, 0])[105:134, 90]) 
    58 ; 11/ Island 
    59 ;   mu[9] = - total( (hv*e1v*vmaskr[*, *, 0])[117:130, 124]) 
    60    mu[11] = - total( (hv*e1v*vmaskr[*, *, 0])[132:143, 120]) 
    61 ; 12/ Spitzberg 
    62 ;   mu[10] = - total( (hv*e1v*vmaskr[*, *, 0])[110:120, 137]) 
    63    mu[12] = - total( (hv*e1v*vmaskr[*, *, 0])[140:160, 135]) 
    64 ; 13/ Japan 
    65    mu[13] = total( (hv*e1v*vmaskr[*, *, 0])[25:131, 100]) 
    66 end 
    67 ;+ 
    68 ; 
    69 ; Lecture du fichier contenant les fonctions de courant associe a 
    70 ; chaque ile (pour ORCA) 
    71 ; 
    72 ; @history 
    73 ; Creation : printemps 1998 / G. Roullet 
    74 ; modification : juin 99 
    75 ; 
    76 ;- 
    77 PRO diag_bsf_lec, bsf 
    78 ; 
    79 @common 
    80 @com_eg 
    81    ;; 
    82    jpisl=14 
    83    bsf = fltarr(jpi, jpj, jpisl) 
    84    z2d = fltarr(jpi, jpj) 
    85    fichier = 'islands.nc' 
    86    FOR k = 1, jpisl DO BEGIN  
    87       ile = 'island'+strtrim(string(k), 2) 
    88       z2d = nc_get(fichier, ile, NCDF_DB =hom_idl+'grids/') 
    89       bsf[*, *, k-1] = shift(z2d, 0, 0) 
    90    ENDFOR  
    91 end 
    92  
    93 ;+ 
    9424; 
    9525;- 
    … …  
    10232@common 
    10333COMMON bsf_lec, bsfislands 
    104  
     34; 
    10535   jpisl=14 
    10636   IF keyword_set(flux) EQ 0 THEN BEGIN  

  • trunk/tools/bsf/diag_bsf_lec.pro

    r170 r173  
    11;+ 
    22; 
    3 ; Calcul de la Fonction de Courant Barotrope  
    4 ; a partir du champ de vitesse (u,v) 
     3; Lecture du fichier contenant les fonctions de courant associe a 
     4; chaque ile (pour ORCA) 
     5; 
     6; @param bsf 
    57; 
    68; @uses 
    … …  
    911; 
    1012; @history 
    11 ; creation : printemps 98 G. Roullet 
    1213; 
    13 ; @version 
    14 ; $Id$ 
     14; - fplod 20091118T075330Z  
     15; 
     16;   * extract from diag_bsf.pro 
     17; 
     18; - modification : juin 99 
     19; 
     20; - Creation : printemps 1998 / G. Roullet 
    1521; 
    1622;- 
    17 PRO diag_bsf_mu, hu, hv, mu 
     23PRO diag_bsf_lec, bsf 
    1824; 
    1925  compile_opt idl2, strictarrsubs 
    … …  
    2228@com_eg 
    2329; 
    24 ;; calcul adapte a la grille ORCA 
    25 ;; 
    26    mu = fltarr(14) 
    27 ; 0/ Antartic 
    28 ;   mu[0] = total( (hu*e2u*umaskr[*, *, 0])[149, 5:35]) 
    29    mu[0] = total( (hu*e2u*umaskr[*, *, 0])[151, 14:48]) 
    30 ; 1/ America 
    31 ;   mu[1] = -total( (hv*e1v*vmaskr[*, *, 0])[110:150, 60]) 
    32    mu[1] = -total( (hv*e1v*vmaskr[*, *, 0])[110:150, 60]) 
    33 ; 3/ Australia 
    34 ;   mu[3] = total( (hu*e2u*umaskr[*, *, 0])[17, 40:49]) 
    35    mu[3] = total( (hu*e2u*umaskr[*, *, 0])[20, 50:61]) 
    36 ; 2/ New Zealand 
    37 ;   mu[2] = mu[3]+ total( (hv*e1v*vmaskr[*, *, 0])[30:40, 25]) 
    38    mu[2] = mu[3]+ total( (hv*e1v*vmaskr[*, *, 0])[35:49, 44]) 
    39 ; 4/ Madagascar 
    40 ;   mu[4] = total( (hv*e1v*vmaskr[*, *, 0])[150:159, 40]) 
    41    mu[4] = total( (hv*e1v*vmaskr[*, *, 0])[160:164, 56]) 
    42 ; 5/ New Guinea 
    43 ;   mu[5] = - total( (hu*e2u*umaskr[*, *, 0])[16, 50:55]) 
    44    mu[5] = total( (hv*e2v*vmaskr[*, *, 0])[23:33, 61]) 
    45 ; 6/ Celebes 
    46    mu[6] = total( (hv*e1v*vmaskr[*, *, 0])[13:22, 64]) 
    47 ; 7/ Borneo 
    48 ;   mu[6] = total( (hv*e1v*vmaskr[*, *, 0])[7:12, 60]) 
    49    mu[7] = total( (hv*e1v*vmaskr[*, *, 0])[14:17, 67]) 
    50 ; 8/ Philippines 
    51 ;   mu[7] = total( (hv*e1v*vmaskr[*, *, 0])[10:16, 90]) 
    52    mu[8] = total( (hv*e1v*vmaskr[*, *, 0])[14:22, 86]) 
    53 ; 9/ Tasmania 
    54    mu[9] = mu[3]+ total( (hu*e1u*umaskr[*, *, 0])[34, 41:45]) 
    55 ; 10/ Cuba 
    56 ;   mu[8] = - total( (hv*e1v*vmaskr[*, *, 0])[101:127, 86]) 
    57    mu[10] = - total( (hv*e1v*vmaskr[*, *, 0])[105:134, 90]) 
    58 ; 11/ Island 
    59 ;   mu[9] = - total( (hv*e1v*vmaskr[*, *, 0])[117:130, 124]) 
    60    mu[11] = - total( (hv*e1v*vmaskr[*, *, 0])[132:143, 120]) 
    61 ; 12/ Spitzberg 
    62 ;   mu[10] = - total( (hv*e1v*vmaskr[*, *, 0])[110:120, 137]) 
    63    mu[12] = - total( (hv*e1v*vmaskr[*, *, 0])[140:160, 135]) 
    64 ; 13/ Japan 
    65    mu[13] = total( (hv*e1v*vmaskr[*, *, 0])[25:131, 100]) 
    66 end 
    67 ;+ 
    68 ; 
    69 ; Lecture du fichier contenant les fonctions de courant associe a 
    70 ; chaque ile (pour ORCA) 
    71 ; 
    72 ; @history 
    73 ; Creation : printemps 1998 / G. Roullet 
    74 ; modification : juin 99 
    75 ; 
    76 ;- 
    77 PRO diag_bsf_lec, bsf 
    78 ; 
    79 @common 
    80 @com_eg 
    81    ;; 
    8230   jpisl=14 
    8331   bsf = fltarr(jpi, jpj, jpisl) 
    … …  
    9038   ENDFOR  
    9139end 
    92  
    93 ;+ 
    94 ; 
    95 ;- 
    96 FUNCTION diag_bsf, u, v $ 
    97          , FLUX = flux $ 
    98          , _EXTRA=extra 
    99 ; 
    100   compile_opt idl2, strictarrsubs 
    101 ; 
    102 @common 
    103 COMMON bsf_lec, bsfislands 
    104  
    105    jpisl=14 
    106    IF keyword_set(flux) EQ 0 THEN BEGIN  
    107 ; 
    108 ;  integration sur la verticale 
    109 ; 
    110       hu = moyz(u, /int) 
    111       hv = moyz(v, /int) 
    112    ENDIF ELSE BEGIN  
    113 ; 
    114 ;  on suppose que u et v sont deja integres verticalement 
    115 ; 
    116       hu = u 
    117       hv = v 
    118    ENDELSE  
    119     
    120  print, '     Inversion...' 
    121    bsf = inversion( curl2(hu, hv), /f, _EXTRA=extra) 
    122 ; 
    123 ;  calcul des potentiels entre chaque ile 
    124 ; 
    125  print, '     Potential...' 
    126    diag_bsf_mu, hu, hv, mu 
    127 ; 
    128 ;  lecture des fonctions courants unite autour de chaque ile 
    129 ; 
    130    IF n_elements(bsfislands) EQ 0 THEN diag_bsf_lec, bsfislands 
    131     
    132    FOR k = 0, jpisl-1 DO bsf = bsf+mu[k]*bsfislands[*, *, k] 
    133  
    134    return, bsf 
    135  
    136 END 

  • trunk/tools/bsf/diag_bsf_mu.pro

    r170 r173  
    11;+ 
    22; 
    3 ; Calcul de la Fonction de Courant Barotrope  
     3; Calcul de la Fonction de Courant Barotrope 
    44; a partir du champ de vitesse (u,v) 
     5; 
     6; @param hu 
     7; 
     8; @param hv 
     9; 
     10; @param mu 
    511; 
    612; @uses 
    … …  
    915; 
    1016; @history 
    11 ; creation : printemps 98 G. Roullet 
     17; 
     18; - fplod 20091118T075330Z  
     19; 
     20;   * extract from diag_bsf.pro 
     21; 
     22; - creation : printemps 98 G. Roullet 
    1223; 
    1324; @version 
    … …  
    6576   mu[13] = total( (hv*e1v*vmaskr[*, *, 0])[25:131, 100]) 
    6677end 
    67 ;+ 
    68 ; 
    69 ; Lecture du fichier contenant les fonctions de courant associe a 
    70 ; chaque ile (pour ORCA) 
    71 ; 
    72 ; @history 
    73 ; Creation : printemps 1998 / G. Roullet 
    74 ; modification : juin 99 
    75 ; 
    76 ;- 
    77 PRO diag_bsf_lec, bsf 
    78 ; 
    79 @common 
    80 @com_eg 
    81    ;; 
    82    jpisl=14 
    83    bsf = fltarr(jpi, jpj, jpisl) 
    84    z2d = fltarr(jpi, jpj) 
    85    fichier = 'islands.nc' 
    86    FOR k = 1, jpisl DO BEGIN  
    87       ile = 'island'+strtrim(string(k), 2) 
    88       z2d = nc_get(fichier, ile, NCDF_DB =hom_idl+'grids/') 
    89       bsf[*, *, k-1] = shift(z2d, 0, 0) 
    90    ENDFOR  
    91 end 
    92  
    93 ;+ 
    94 ; 
    95 ;- 
    96 FUNCTION diag_bsf, u, v $ 
    97          , FLUX = flux $ 
    98          , _EXTRA=extra 
    99 ; 
    100   compile_opt idl2, strictarrsubs 
    101 ; 
    102 @common 
    103 COMMON bsf_lec, bsfislands 
    104  
    105    jpisl=14 
    106    IF keyword_set(flux) EQ 0 THEN BEGIN  
    107 ; 
    108 ;  integration sur la verticale 
    109 ; 
    110       hu = moyz(u, /int) 
    111       hv = moyz(v, /int) 
    112    ENDIF ELSE BEGIN  
    113 ; 
    114 ;  on suppose que u et v sont deja integres verticalement 
    115 ; 
    116       hu = u 
    117       hv = v 
    118    ENDELSE  
    119     
    120  print, '     Inversion...' 
    121    bsf = inversion( curl2(hu, hv), /f, _EXTRA=extra) 
    122 ; 
    123 ;  calcul des potentiels entre chaque ile 
    124 ; 
    125  print, '     Potential...' 
    126    diag_bsf_mu, hu, hv, mu 
    127 ; 
    128 ;  lecture des fonctions courants unite autour de chaque ile 
    129 ; 
    130    IF n_elements(bsfislands) EQ 0 THEN diag_bsf_lec, bsfislands 
    131     
    132    FOR k = 0, jpisl-1 DO bsf = bsf+mu[k]*bsfislands[*, *, k] 
    133  
    134    return, bsf 
    135  
    136 END 

  • trunk/tools/msf_bin/bin_velocity.pro

    r170 r173  
    11;+ 
    22; 
    3 ; Calcule la fonction de courant meridienne definie en coordonnees 
    4 ; isopycnes 
    5 ; i.e. calcule MSF(lat,sigma) au lieu de MSF(lat,z) 
     3; @param dens 
    64; 
    7 ; ATTENTION ! u_s et v_s contiennent le transport par 
    8 ; unité de longueur, i.e il s'agit de u_s=u.dz intégré 
    9 ; sur la couche d'isopycne considérée et non pas de la 
    10 ; vitesse moyenne dans la couche 
     5; @param sigma_values 
    116; 
    12 ; @param DENSITY {in}{required} 
     7; @param x1 
    138; 
    14 ; @param SIGMA_VALUES {in}{required} 
    15 ; valeurs de bining 
     9; @keyword UU 
    1610; 
    17 ; @param U {in}{required} 
     11; @keyword VV 
    1812; 
    19 ; @param V {in}{required} 
     13; @uses 
     14; <pro>common</pro> 
    2015; 
    2116; @history 
    22 ; 2001/06/05 G. Roullet <roullet\@univ-brest.fr> 
     17; - fplod 20091118T084527Z 
     18; 
     19;   * extract from msf_sigma.pro 
     20; 
     21; - 2001/06/05 G. Roullet <roullet\@univ-brest.fr> 
    2322; 
    2423; @version 
    … …  
    2625; 
    2726;- 
    28 FUNCTION msf_sigma, density, sigma_values, u, v $ 
    29          , MSK=msk 
    30 ; 
    31   compile_opt idl2, strictarrsubs 
    32 ; 
    33       print, 'Bining U ...' 
    34       u_s = bin_velocity(density, sigma_values, u, /uu) 
    35  
    36       print, 'Bining V ...' 
    37       v_s = bin_velocity(density, sigma_values, v, /vv) 
    38  
    39       print, 'Computing MSF ...' 
    40       compute_msf, u_s, v_s, msf, MSK = msk 
    41  
    42       return, msf 
    43 END   
    44  
    45 ;+ 
    46 ; 
    47 ; @uses 
    48 ; <pro>common</pro> 
    49 ; 
    50 ;- 
    51 FUNCTION bin_velocity, dens, sigma_values, x1 $  
    52          , UU = uu $ 
    53          , VV = vv 
     27FUNCTION bin_velocity, dens, sigma_values, x1 $ 
     28         , UU=uu $ 
     29         , VV=vv 
    5430; 
    5531  compile_opt idl2, strictarrsubs 
    5632; 
    5733@common 
    58     
     34; 
    5935   IF keyword_set(uu) THEN BEGIN 
    6036      density = boundperio(umask*(dens+shift(dens, -1, 0, 0))/2, /uu) 
    6137      xmask = umask 
    62    ENDIF  
     38   ENDIF 
    6339 
    6440   IF keyword_set(vv) THEN BEGIN 
    6541      density = boundperio(vmask*(dens+shift(dens, 0, -1, 0))/2, /vv) 
    6642      xmask = vmask 
    67    ENDIF  
     43   ENDIF 
    6844 
    6945      s_s = sigma_values 
    7046      N_s = n_elements(s_s) 
    7147      N_z = jpk 
    72        
     48 
    7349;; 
    7450;; Definition des variables 
    … …  
    9672   x1_content = fltarr(jpi, jpj, jpk) 
    9773   x1_content[*, *, 0] = e3t[0] * x1[*, *, 0] 
    98    FOR k = 1, jpk-1 DO x1_content[*, *, k] = x1_content[*, *, k-1] $  
     74   FOR k = 1, jpk-1 DO x1_content[*, *, k] = x1_content[*, *, k-1] $ 
    9975                        + e3t[k] * x1[*, *, k] 
    10076 
    … …  
    10278;  Boucle sur le domaine 2D 
    10379; 
    104    FOR i = 0, (jpi-1) DO BEGIN  
    105       FOR j = 0, (jpj-1) DO BEGIN  
     80   FOR i = 0, (jpi-1) DO BEGIN 
     81      FOR j = 0, (jpj-1) DO BEGIN 
    10682; 
    10783;  Indices des points T dans l ocean 
    … …  
    11490 
    11591         IF (i_ocean[0] NE -1) THEN BEGIN ; on n'entre que si il y a des points ocean 
    116   
     92 
    11793            s_z[*] = density[i, j, *] 
    11894            c1_z[*] = x1_content[i, j, *] 
    … …  
    12197; pathologiques en attendant d'ecrire une vraie commande d'extraction 
    12298; de progils stritement croissant. Il s'agit donc d'un test adhoc. 
    123             IF s_z[0] LT s_z[i_ocean[n_elements(i_ocean)-1]] THEN BEGIN  
     99            IF s_z[0] LT s_z[i_ocean[n_elements(i_ocean)-1]] THEN BEGIN 
    124100 
    125101;------------------------------------------------------------------------ 
    … …  
    137113            z_s[N_s) = z_zw(i_bottom] 
    138114            c1_s[N_s] = x1_content[i, j, jpk-1] 
    139              
     115 
    140116; extraction d'un sous profil strictement croissant 
    141117            mini = min(s_z[i_ocean]) 
    142             maxi = max(s_z[i_ocean])             
     118            maxi = max(s_z[i_ocean]) 
    143119 
    144120            i_min = where(s_z[i_ocean] EQ mini) 
    … …  
    151127;            IF i_max GE jpk-1 THEN print, i, j, i_max 
    152128;            IF i_min LE 1 THEN print, i, j, i_min 
    153              
    154 ; Si la valeur du niveau (s_s) est plus faible que la densite de surface,  
     129 
     130; Si la valeur du niveau (s_s) est plus faible que la densite de surface, 
    155131; l isopycne est mise en surface (z_s=0) 
    156132; 
    157133            ind = where(s_s LT s_z[i_min]) 
    158             IF ind[0] NE -1 THEN BEGIN  
     134            IF ind[0] NE -1 THEN BEGIN 
    159135               z_s[ind] = 0 
    160136               c1_s[ind] = 0 
    161             ENDIF  
     137            ENDIF 
    162138 
    163 ; Si la valeur du niveau (s_s) est plus elevee que la densite du fond,  
     139; Si la valeur du niveau (s_s) est plus elevee que la densite du fond, 
    164140; l isopycne est mise au fond (z_s=z_zw(i_bottom)) 
    165141; 
    166142            ind = where(s_s GT s_z[i_max]) 
    167             IF ind[0] NE -1 THEN BEGIN  
     143            IF ind[0] NE -1 THEN BEGIN 
    168144               z_s[ind] = z_s[N_s] 
    169145               c1_s[ind] = c1_s[N_s] 
    170             ENDIF  
     146            ENDIF 
    171147 
    172148            ind = where( (s_s GE s_z[i_min]) AND (s_s LT s_z[i_max]) ) 
    173             IF ind[0] NE -1 THEN BEGIN  
    174                 
     149            IF ind[0] NE -1 THEN BEGIN 
     150 
    175151               i_profil = i_ocean[i_min:i_max] 
    176152 
    177153               z_s[ind] = interpol(z_zt[i_profil], s_z[i_profil], s_s[ind]) 
    178                 
     154 
    179155; 
    180156; j'utilise la fonction spline pour interpoler le contenu 
    … …  
    191167               x1_s[ind] = (c1_s[ind+1]-c1_s[ind]);/(z_s[ind+1]-z_s[ind]) 
    192168 
    193             ENDIF  
    194          ENDIF  
     169            ENDIF 
     170         ENDIF 
    195171         ENDIF 
    196172 
    197173         x1_bin[i, j, *] = x1_s 
    198           
    199       ENDFOR  
    200    ENDFOR  
    201     
     174 
     175      ENDFOR 
     176   ENDFOR 
     177 
    202178   return, x1_bin 
    203179 
    204 END  
    205 ;+ 
    206 ; 
    207 ; @uses 
    208 ; <pro>common</pro> 
    209 ; 
    210 ;- 
    211 PRO compute_msf, u, v, msf $ 
    212     , MSK = msk 
    213 ; 
    214   compile_opt idl2, strictarrsubs 
    215 ; 
    216 ; 
    217 @common 
    218 COMMON flx, indxu, indxv, indyu, indyv, zxuleng, zxvleng, zyuleng, zyvleng 
    219     
    220    IF n_elements(indxu) EQ 0 THEN iniflx 
    221  
    222  
    223 ; 
    224 ;  jpl bandes de latitudes 
    225 ; 
    226    jpk2 = (size(u))(3) ; on remplace jpk par jpk2 = nombre de couches 
    227  
    228  
    229    fm = fltarr(jpl, jpk2) 
    230    msf = fltarr(jpl, jpk2) 
    231    msfmsk = fltarr(jpl, jpk2) 
    232  
    233    vert = replicate(1, jpk2) 
    234  
    235    zu = u 
    236    zv = v 
    237 ; 
    238 ;  Masque de u et v pour un calcul de msf sur le sous domaine msk(jpi,jpj) 
    239 ; 
    240    IF n_elements(msk) NE 0 THEN BEGIN  
    241       zumask = boundperio( (msk +shift(msk, -1, 0) ) < 1 ) 
    242       zvmask = boundperio( (msk +shift(msk, 0, -1) ) < 1 ) 
    243       zumask = zumask[*]#vert 
    244       zvmask = zvmask[*]#vert 
    245       zu = zu*zumask 
    246       zv = zv*zvmask 
    247    ENDIF  
    248 ; 
    249 ;  Ecriture "tricky" optimisee... si vous avez compris comment 
    250 ;  on calculait un flux a partir de champ 2D, vous avez fait le 
    251 ;  plus dur : ici on generalise a 3D d''ou l''utilisation de #vert 
    252 ;  
    253 ; 
    254 ;  indxu pointe dans le tableau [jpi,jpj], il faut donc l''etendre 
    255 ;  pour pointer dans le tableau [jpi,jpj,jpk2] et ne pas oublier 
    256 ;  le decalage : le niveau k a ses indices decales de k*jpi*jpj par 
    257 ;  rapport au niveau 0 
    258 ; 
    259    offset = (fltarr(180, jpl)+long(jpi)*jpj)(*)#findgen(jpk2) 
    260    indu = indxu[*]#vert+offset 
    261    indv = indxv[*]#vert+offset    
    262 ; 
    263 ;  extension sur la verticale 
    264 ; 
    265    zuleng = zxuleng[*]#vert 
    266    zvleng = zxvleng[*]#vert 
    267 ; 
    268 ;  calcul du flux 
    269 ; 
    270    z = reform(zu[indu]*zuleng+zv[indv]*zvleng, 180, jpl, jpk2) 
    271 ; 
    272 ;  integration zonale du flux ! 
    273 ; 
    274    fm = total(z, 1)  
    275 ; 
    276 ;  calcul de la msf en integrant depuis le fond 
    277 ; 
    278 ;        Remarque : on ne multiplie pas par e3t car les champs 
    279 ;        d'entree sont deja integrés verticalement (sur l'epaisseur de 
    280 ;        couche). 
    281 ; 
    282    FOR k = jpk2-2, 0, -1 DO msf[*, k] = msf[*, k+1]-fm[*, k]; *e3t[k] 
    283  
    284180END 

  • trunk/tools/msf_bin/compute_msf.pro

    r170 r173  
    11;+ 
    22; 
    3 ; Calcule la fonction de courant meridienne definie en coordonnees 
    4 ; isopycnes 
    5 ; i.e. calcule MSF(lat,sigma) au lieu de MSF(lat,z) 
     3; @param U 
    64; 
    7 ; ATTENTION ! u_s et v_s contiennent le transport par 
    8 ; unité de longueur, i.e il s'agit de u_s=u.dz intégré 
    9 ; sur la couche d'isopycne considérée et non pas de la 
    10 ; vitesse moyenne dans la couche 
     5; @param V 
    116; 
    12 ; @param DENSITY {in}{required} 
     7; @param MSF 
    138; 
    14 ; @param SIGMA_VALUES {in}{required} 
    15 ; valeurs de bining 
     9; @keyword MSK 
    1610; 
    17 ; @param U {in}{required} 
    18 ; 
    19 ; @param V {in}{required} 
     11; @uses 
     12; <pro>common</pro> 
    2013; 
    2114; @history 
     15; - fplod 20091118T084809Z 
     16; 
     17;   * extract from msf_sigma.pro 
     18; 
    2219; 2001/06/05 G. Roullet <roullet\@univ-brest.fr> 
    2320; 
    2421; @version 
    2522; $Id$ 
    26 ; 
    27 ;- 
    28 FUNCTION msf_sigma, density, sigma_values, u, v $ 
    29          , MSK=msk 
    30 ; 
    31   compile_opt idl2, strictarrsubs 
    32 ; 
    33       print, 'Bining U ...' 
    34       u_s = bin_velocity(density, sigma_values, u, /uu) 
    35  
    36       print, 'Bining V ...' 
    37       v_s = bin_velocity(density, sigma_values, v, /vv) 
    38  
    39       print, 'Computing MSF ...' 
    40       compute_msf, u_s, v_s, msf, MSK = msk 
    41  
    42       return, msf 
    43 END   
    44  
    45 ;+ 
    46 ; 
    47 ; @uses 
    48 ; <pro>common</pro> 
    49 ; 
    50 ;- 
    51 FUNCTION bin_velocity, dens, sigma_values, x1 $  
    52          , UU = uu $ 
    53          , VV = vv 
    54 ; 
    55   compile_opt idl2, strictarrsubs 
    56 ; 
    57 @common 
    58     
    59    IF keyword_set(uu) THEN BEGIN 
    60       density = boundperio(umask*(dens+shift(dens, -1, 0, 0))/2, /uu) 
    61       xmask = umask 
    62    ENDIF  
    63  
    64    IF keyword_set(vv) THEN BEGIN 
    65       density = boundperio(vmask*(dens+shift(dens, 0, -1, 0))/2, /vv) 
    66       xmask = vmask 
    67    ENDIF  
    68  
    69       s_s = sigma_values 
    70       N_s = n_elements(s_s) 
    71       N_z = jpk 
    72        
    73 ;; 
    74 ;; Definition des variables 
    75 ;; 
    76 ;   Profils selon les niveaux du modele (suffixe _z) 
    77 ; 
    78    c1_z = fltarr(N_z)    ; profil du contenu vertical de x 
    79    s_z = fltarr(N_z)     ; profil de la densite 
    80    z_zt = zdept[1:jpk]  ; profondeur au point T  (k=0 -> 5m) 
    81    z_zw = zdepw[1:jpk]  ;                     W  (k=0 -> 10m) 
    82  
    83 ;   Profils selon les couches de densite (suffixe _s) 
    84 ; 
    85    z_s = fltarr(N_s+1) 
    86    c1_s = fltarr(N_s+1) 
    87    x1_s = fltarr(N_s+1) 
    88  
    89 ; Tableau de sortie 
    90 ; 
    91    x1_bin = fltarr(jpi, jpj, N_s+1) 
    92  
    93 ; 
    94 ;  Calcul du contenu de x sur la verticale (permet d assurer la conservation integrale) 
    95 ; 
    96    x1_content = fltarr(jpi, jpj, jpk) 
    97    x1_content[*, *, 0] = e3t[0] * x1[*, *, 0] 
    98    FOR k = 1, jpk-1 DO x1_content[*, *, k] = x1_content[*, *, k-1] $  
    99                         + e3t[k] * x1[*, *, k] 
    100  
    101 ; 
    102 ;  Boucle sur le domaine 2D 
    103 ; 
    104    FOR i = 0, (jpi-1) DO BEGIN  
    105       FOR j = 0, (jpj-1) DO BEGIN  
    106 ; 
    107 ;  Indices des points T dans l ocean 
    108 ; 
    109          i_ocean = where(xmask[i, j, *] EQ 1) 
    110  
    111          z_s = z_s*0. 
    112          c1_s = c1_s*0. 
    113          x1_s = x1_s*0. 
    114  
    115          IF (i_ocean[0] NE -1) THEN BEGIN ; on n'entre que si il y a des points ocean 
    116   
    117             s_z[*] = density[i, j, *] 
    118             c1_z[*] = x1_content[i, j, *] 
    119  
    120 ; critere supplementaire a imposer sur le profil pour eviter les cas 
    121 ; pathologiques en attendant d'ecrire une vraie commande d'extraction 
    122 ; de progils stritement croissant. Il s'agit donc d'un test adhoc. 
    123             IF s_z[0] LT s_z[i_ocean[n_elements(i_ocean)-1]] THEN BEGIN  
    124  
    125 ;------------------------------------------------------------------------ 
    126 ; controle si le profil est bien strictement croissant (pour l'instant 
    127 ; inutilisé (avis aux amateurs) 
    128 ; 
    129 ;            ds = (shift(s_z, -1)-s_z)[i_ocean[0:n_elements(i_ocean)-2]] 
    130 ;            ds[0] = 0 
    131 ;            ind = where(ds LT 0.) 
    132 ;            croissant = ind[0] EQ -1 
    133 ;------------------------------------------------------------------------ 
    134  
    135             i_bottom = i_ocean[n_elements(i_ocean)-1] 
    136  
    137             z_s[N_s) = z_zw(i_bottom] 
    138             c1_s[N_s] = x1_content[i, j, jpk-1] 
    139              
    140 ; extraction d'un sous profil strictement croissant 
    141             mini = min(s_z[i_ocean]) 
    142             maxi = max(s_z[i_ocean])             
    143  
    144             i_min = where(s_z[i_ocean] EQ mini) 
    145             i_max = where(s_z[i_ocean] EQ maxi) 
    146 ;   on prend le plus grand des indices min 
    147 ;         et le plus petit des indices max 
    148             i_min = i_min(n_elements(i_min)-1) 
    149             i_max = i_max[0] 
    150  
    151 ;            IF i_max GE jpk-1 THEN print, i, j, i_max 
    152 ;            IF i_min LE 1 THEN print, i, j, i_min 
    153              
    154 ; Si la valeur du niveau (s_s) est plus faible que la densite de surface,  
    155 ; l isopycne est mise en surface (z_s=0) 
    156 ; 
    157             ind = where(s_s LT s_z[i_min]) 
    158             IF ind[0] NE -1 THEN BEGIN  
    159                z_s[ind] = 0 
    160                c1_s[ind] = 0 
    161             ENDIF  
    162  
    163 ; Si la valeur du niveau (s_s) est plus elevee que la densite du fond,  
    164 ; l isopycne est mise au fond (z_s=z_zw(i_bottom)) 
    165 ; 
    166             ind = where(s_s GT s_z[i_max]) 
    167             IF ind[0] NE -1 THEN BEGIN  
    168                z_s[ind] = z_s[N_s] 
    169                c1_s[ind] = c1_s[N_s] 
    170             ENDIF  
    171  
    172             ind = where( (s_s GE s_z[i_min]) AND (s_s LT s_z[i_max]) ) 
    173             IF ind[0] NE -1 THEN BEGIN  
    174                 
    175                i_profil = i_ocean[i_min:i_max] 
    176  
    177                z_s[ind] = interpol(z_zt[i_profil], s_z[i_profil], s_s[ind]) 
    178                 
    179 ; 
    180 ; j'utilise la fonction spline pour interpoler le contenu 
    181 ; 
    182                c1_s[ind] = spline(z_zw[i_profil], c1_z[i_profil], z_s[ind], 1) 
    183 ; 
    184 ; l'interpolation lineaire marche aussi. Elle donne des resultats differents. Elle 
    185 ; me semble moins propre. Je la donne en remarque. 
    186 ; 
    187 ;               c_s[ind+1] = interpol(c_z[i_profil], z_zw[i_profil], z_s[ind+1]) 
    188 ; 
    189 ; Remarque : on ne divise pas par l'epaisseur de la couche 
    190 ; 
    191                x1_s[ind] = (c1_s[ind+1]-c1_s[ind]);/(z_s[ind+1]-z_s[ind]) 
    192  
    193             ENDIF  
    194          ENDIF  
    195          ENDIF 
    196  
    197          x1_bin[i, j, *] = x1_s 
    198           
    199       ENDFOR  
    200    ENDFOR  
    201     
    202    return, x1_bin 
    203  
    204 END  
    205 ;+ 
    206 ; 
    207 ; @uses 
    208 ; <pro>common</pro> 
    20923; 
    21024;- 
    … …  
    21731@common 
    21832COMMON flx, indxu, indxv, indyu, indyv, zxuleng, zxvleng, zyuleng, zyvleng 
    219     
     33 
    22034   IF n_elements(indxu) EQ 0 THEN iniflx 
    22135 
    … …  
    22438;  jpl bandes de latitudes 
    22539; 
    226    jpk2 = (size(u))(3) ; on remplace jpk par jpk2 = nombre de couches 
     40   jpk2 = (size(u))[3] ; on remplace jpk par jpk2 = nombre de couches 
    22741 
    22842 
    … …  
    23852;  Masque de u et v pour un calcul de msf sur le sous domaine msk(jpi,jpj) 
    23953; 
    240    IF n_elements(msk) NE 0 THEN BEGIN  
     54   IF n_elements(msk) NE 0 THEN BEGIN 
    24155      zumask = boundperio( (msk +shift(msk, -1, 0) ) < 1 ) 
    24256      zvmask = boundperio( (msk +shift(msk, 0, -1) ) < 1 ) 
    … …  
    24559      zu = zu*zumask 
    24660      zv = zv*zvmask 
    247    ENDIF  
     61   ENDIF 
    24862; 
    24963;  Ecriture "tricky" optimisee... si vous avez compris comment 
    25064;  on calculait un flux a partir de champ 2D, vous avez fait le 
    25165;  plus dur : ici on generalise a 3D d''ou l''utilisation de #vert 
    252 ;  
     66; 
    25367; 
    25468;  indxu pointe dans le tableau [jpi,jpj], il faut donc l''etendre 
    … …  
    25771;  rapport au niveau 0 
    25872; 
    259    offset = (fltarr(180, jpl)+long(jpi)*jpj)(*)#findgen(jpk2) 
     73   offset = (fltarr(180, jpl)+long(jpi)*jpj)[*]#findgen(jpk2) 
    26074   indu = indxu[*]#vert+offset 
    261    indv = indxv[*]#vert+offset    
     75   indv = indxv[*]#vert+offset 
    26276; 
    26377;  extension sur la verticale 
    … …  
    27286;  integration zonale du flux ! 
    27387; 
    274    fm = total(z, 1)  
     88   fm = total(z, 1) 
    27589; 
    27690;  calcul de la msf en integrant depuis le fond 

  • trunk/tools/msf_bin/msf_sigma.pro

    r170 r173  
    1010; vitesse moyenne dans la couche 
    1111; 
    12 ; @param DENSITY {in}{required} 
     12; @param density {in}{required} 
    1313; 
    14 ; @param SIGMA_VALUES {in}{required} 
     14; @param sigma_values {in}{required} 
    1515; valeurs de bining 
    1616; 
    17 ; @param U {in}{required} 
     17; @param u {in}{required} 
    1818; 
    19 ; @param V {in}{required} 
     19; @param v {in}{required} 
     20; 
     21; @keyword MASK 
    2022; 
    2123; @history 
    22 ; 2001/06/05 G. Roullet <roullet\@univ-brest.fr> 
     24; - fplod 20091118T085107Z  
     25;    
     26;   * externalize bin_velocity and compute_msf 
     27; 
     28; - 2001/06/05 G. Roullet <roullet\@univ-brest.fr> 
    2329; 
    2430; @version 
    … …  
    4248      return, msf 
    4349END   
    44  
    45 ;+ 
    46 ; 
    47 ; @uses 
    48 ; <pro>common</pro> 
    49 ; 
    50 ;- 
    51 FUNCTION bin_velocity, dens, sigma_values, x1 $  
    52          , UU = uu $ 
    53          , VV = vv 
    54 ; 
    55   compile_opt idl2, strictarrsubs 
    56 ; 
    57 @common 
    58     
    59    IF keyword_set(uu) THEN BEGIN 
    60       density = boundperio(umask*(dens+shift(dens, -1, 0, 0))/2, /uu) 
    61       xmask = umask 
    62    ENDIF  
    63  
    64    IF keyword_set(vv) THEN BEGIN 
    65       density = boundperio(vmask*(dens+shift(dens, 0, -1, 0))/2, /vv) 
    66       xmask = vmask 
    67    ENDIF  
    68  
    69       s_s = sigma_values 
    70       N_s = n_elements(s_s) 
    71       N_z = jpk 
    72        
    73 ;; 
    74 ;; Definition des variables 
    75 ;; 
    76 ;   Profils selon les niveaux du modele (suffixe _z) 
    77 ; 
    78    c1_z = fltarr(N_z)    ; profil du contenu vertical de x 
    79    s_z = fltarr(N_z)     ; profil de la densite 
    80    z_zt = zdept[1:jpk]  ; profondeur au point T  (k=0 -> 5m) 
    81    z_zw = zdepw[1:jpk]  ;                     W  (k=0 -> 10m) 
    82  
    83 ;   Profils selon les couches de densite (suffixe _s) 
    84 ; 
    85    z_s = fltarr(N_s+1) 
    86    c1_s = fltarr(N_s+1) 
    87    x1_s = fltarr(N_s+1) 
    88  
    89 ; Tableau de sortie 
    90 ; 
    91    x1_bin = fltarr(jpi, jpj, N_s+1) 
    92  
    93 ; 
    94 ;  Calcul du contenu de x sur la verticale (permet d assurer la conservation integrale) 
    95 ; 
    96    x1_content = fltarr(jpi, jpj, jpk) 
    97    x1_content[*, *, 0] = e3t[0] * x1[*, *, 0] 
    98    FOR k = 1, jpk-1 DO x1_content[*, *, k] = x1_content[*, *, k-1] $  
    99                         + e3t[k] * x1[*, *, k] 
    100  
    101 ; 
    102 ;  Boucle sur le domaine 2D 
    103 ; 
    104    FOR i = 0, (jpi-1) DO BEGIN  
    105       FOR j = 0, (jpj-1) DO BEGIN  
    106 ; 
    107 ;  Indices des points T dans l ocean 
    108 ; 
    109          i_ocean = where(xmask[i, j, *] EQ 1) 
    110  
    111          z_s = z_s*0. 
    112          c1_s = c1_s*0. 
    113          x1_s = x1_s*0. 
    114  
    115          IF (i_ocean[0] NE -1) THEN BEGIN ; on n'entre que si il y a des points ocean 
    116   
    117             s_z[*] = density[i, j, *] 
    118             c1_z[*] = x1_content[i, j, *] 
    119  
    120 ; critere supplementaire a imposer sur le profil pour eviter les cas 
    121 ; pathologiques en attendant d'ecrire une vraie commande d'extraction 
    122 ; de progils stritement croissant. Il s'agit donc d'un test adhoc. 
    123             IF s_z[0] LT s_z[i_ocean[n_elements(i_ocean)-1]] THEN BEGIN  
    124  
    125 ;------------------------------------------------------------------------ 
    126 ; controle si le profil est bien strictement croissant (pour l'instant 
    127 ; inutilisé (avis aux amateurs) 
    128 ; 
    129 ;            ds = (shift(s_z, -1)-s_z)[i_ocean[0:n_elements(i_ocean)-2]] 
    130 ;            ds[0] = 0 
    131 ;            ind = where(ds LT 0.) 
    132 ;            croissant = ind[0] EQ -1 
    133 ;------------------------------------------------------------------------ 
    134  
    135             i_bottom = i_ocean[n_elements(i_ocean)-1] 
    136  
    137             z_s[N_s) = z_zw(i_bottom] 
    138             c1_s[N_s] = x1_content[i, j, jpk-1] 
    139              
    140 ; extraction d'un sous profil strictement croissant 
    141             mini = min(s_z[i_ocean]) 
    142             maxi = max(s_z[i_ocean])             
    143  
    144             i_min = where(s_z[i_ocean] EQ mini) 
    145             i_max = where(s_z[i_ocean] EQ maxi) 
    146 ;   on prend le plus grand des indices min 
    147 ;         et le plus petit des indices max 
    148             i_min = i_min(n_elements(i_min)-1) 
    149             i_max = i_max[0] 
    150  
    151 ;            IF i_max GE jpk-1 THEN print, i, j, i_max 
    152 ;            IF i_min LE 1 THEN print, i, j, i_min 
    153              
    154 ; Si la valeur du niveau (s_s) est plus faible que la densite de surface,  
    155 ; l isopycne est mise en surface (z_s=0) 
    156 ; 
    157             ind = where(s_s LT s_z[i_min]) 
    158             IF ind[0] NE -1 THEN BEGIN  
    159                z_s[ind] = 0 
    160                c1_s[ind] = 0 
    161             ENDIF  
    162  
    163 ; Si la valeur du niveau (s_s) est plus elevee que la densite du fond,  
    164 ; l isopycne est mise au fond (z_s=z_zw(i_bottom)) 
    165 ; 
    166             ind = where(s_s GT s_z[i_max]) 
    167             IF ind[0] NE -1 THEN BEGIN  
    168                z_s[ind] = z_s[N_s] 
    169                c1_s[ind] = c1_s[N_s] 
    170             ENDIF  
    171  
    172             ind = where( (s_s GE s_z[i_min]) AND (s_s LT s_z[i_max]) ) 
    173             IF ind[0] NE -1 THEN BEGIN  
    174                 
    175                i_profil = i_ocean[i_min:i_max] 
    176  
    177                z_s[ind] = interpol(z_zt[i_profil], s_z[i_profil], s_s[ind]) 
    178                 
    179 ; 
    180 ; j'utilise la fonction spline pour interpoler le contenu 
    181 ; 
    182                c1_s[ind] = spline(z_zw[i_profil], c1_z[i_profil], z_s[ind], 1) 
    183 ; 
    184 ; l'interpolation lineaire marche aussi. Elle donne des resultats differents. Elle 
    185 ; me semble moins propre. Je la donne en remarque. 
    186 ; 
    187 ;               c_s[ind+1] = interpol(c_z[i_profil], z_zw[i_profil], z_s[ind+1]) 
    188 ; 
    189 ; Remarque : on ne divise pas par l'epaisseur de la couche 
    190 ; 
    191                x1_s[ind] = (c1_s[ind+1]-c1_s[ind]);/(z_s[ind+1]-z_s[ind]) 
    192  
    193             ENDIF  
    194          ENDIF  
    195          ENDIF 
    196  
    197          x1_bin[i, j, *] = x1_s 
    198           
    199       ENDFOR  
    200    ENDFOR  
    201     
    202    return, x1_bin 
    203  
    204 END  
    205 ;+ 
    206 ; 
    207 ; @uses 
    208 ; <pro>common</pro> 
    209 ; 
    210 ;- 
    211 PRO compute_msf, u, v, msf $ 
    212     , MSK = msk 
    213 ; 
    214   compile_opt idl2, strictarrsubs 
    215 ; 
    216 ; 
    217 @common 
    218 COMMON flx, indxu, indxv, indyu, indyv, zxuleng, zxvleng, zyuleng, zyvleng 
    219     
    220    IF n_elements(indxu) EQ 0 THEN iniflx 
    221  
    222  
    223 ; 
    224 ;  jpl bandes de latitudes 
    225 ; 
    226    jpk2 = (size(u))(3) ; on remplace jpk par jpk2 = nombre de couches 
    227  
    228  
    229    fm = fltarr(jpl, jpk2) 
    230    msf = fltarr(jpl, jpk2) 
    231    msfmsk = fltarr(jpl, jpk2) 
    232  
    233    vert = replicate(1, jpk2) 
    234  
    235    zu = u 
    236    zv = v 
    237 ; 
    238 ;  Masque de u et v pour un calcul de msf sur le sous domaine msk(jpi,jpj) 
    239 ; 
    240    IF n_elements(msk) NE 0 THEN BEGIN  
    241       zumask = boundperio( (msk +shift(msk, -1, 0) ) < 1 ) 
    242       zvmask = boundperio( (msk +shift(msk, 0, -1) ) < 1 ) 
    243       zumask = zumask[*]#vert 
    244       zvmask = zvmask[*]#vert 
    245       zu = zu*zumask 
    246       zv = zv*zvmask 
    247    ENDIF  
    248 ; 
    249 ;  Ecriture "tricky" optimisee... si vous avez compris comment 
    250 ;  on calculait un flux a partir de champ 2D, vous avez fait le 
    251 ;  plus dur : ici on generalise a 3D d''ou l''utilisation de #vert 
    252 ;  
    253 ; 
    254 ;  indxu pointe dans le tableau [jpi,jpj], il faut donc l''etendre 
    255 ;  pour pointer dans le tableau [jpi,jpj,jpk2] et ne pas oublier 
    256 ;  le decalage : le niveau k a ses indices decales de k*jpi*jpj par 
    257 ;  rapport au niveau 0 
    258 ; 
    259    offset = (fltarr(180, jpl)+long(jpi)*jpj)(*)#findgen(jpk2) 
    260    indu = indxu[*]#vert+offset 
    261    indv = indxv[*]#vert+offset    
    262 ; 
    263 ;  extension sur la verticale 
    264 ; 
    265    zuleng = zxuleng[*]#vert 
    266    zvleng = zxvleng[*]#vert 
    267 ; 
    268 ;  calcul du flux 
    269 ; 
    270    z = reform(zu[indu]*zuleng+zv[indv]*zvleng, 180, jpl, jpk2) 
    271 ; 
    272 ;  integration zonale du flux ! 
    273 ; 
    274    fm = total(z, 1)  
    275 ; 
    276 ;  calcul de la msf en integrant depuis le fond 
    277 ; 
    278 ;        Remarque : on ne multiplie pas par e3t car les champs 
    279 ;        d'entree sont deja integrés verticalement (sur l'epaisseur de 
    280 ;        couche). 
    281 ; 
    282    FOR k = jpk2-2, 0, -1 DO msf[*, k] = msf[*, k+1]-fm[*, k]; *e3t[k] 
    283  
    284 END 

  • trunk/tools/nadw.pro

    r170 r173  
    3838 
    3939end 
    40  
Note: See TracChangeset for help on using the changeset viewer.