;
; mask zonal mean
;
PRO mask_z, fld, cmd, boite_pltz, dimplot, legz

;
; mask pltz data according to bathy or box
;

@common
@com_eg
   legz = ''

   CASE mesh_type OF
      'oce': BEGIN 

         IF splot EQ 1 THEN vertical_domain = [sig_min, sig_max] ELSE vertical_domain = [0, depth_z]

         CASE data_domain OF
            'equator': BEGIN
               tmask[*, *, *] = 1
               IF dimplot EQ 1 THEN boite_pltz = box_h
               IF dimplot EQ 2 THEN boite_pltz = [box_h, vertical_domain]
            END 
            'zonal': BEGIN
                tmask[*, *, *] = 1
                zbox = def_box(cmd.plt, dimplot, legz, time_stride)
                boite_pltz = zbox
                IF dimplot EQ 2 THEN BEGIN
                   IF n_elements(boite_pltz) EQ 4 THEN boite_pltz = [zbox, vertical_domain]
                ENDIF 
;               IF dimplot EQ 1 THEN boite_pltz = box_h
;               IF dimplot EQ 2 THEN boite_pltz = [box_h,
;               vertical_domain]
                idxmsk=where(fld eq valmask)
                IF idxmsk(0) ne -1 THEN BEGIN 
                    print, '     mask_z: building mask from fld valmask'
                    IF (size(fld))(0) EQ 3 THEN BEGIN
                        zmsk = fld(*,*,1) LT valmask/10.
                    ENDIF ELSE BEGIN 
                        zmsk = fld LT valmask/10.
                    ENDELSE 
                    tmask=reform(replicate(jpi,1)#zmsk[*],jpi,jpj,jpk,/overwrite)
                ENDIF
            END  
            ELSE: BEGIN
               ; ORCA geographical domain 
               ; mask with medit
               IF strlen(cmd.plt) EQ 2 THEN BEGIN
                  IF cmd.var NE 'vozonbsf' THEN BEGIN
                     IF (dimplot EQ 1 AND (size(fld))[0] EQ 2) $
                      OR (dimplot EQ 2 AND (size(fld))[0] EQ 3) THEN BEGIN 
                        maskzm = 1-maskread('mediteranean', 'orca', /D3)
                        print, '     [masking out mediteranean]'
                        print, ' '
                        tmask = tmask*maskzm
                     ENDIF 
                     IF dimplot EQ 1 THEN boite_pltz = box_h
                     IF dimplot EQ 2 THEN BEGIN
                        CASE strmid(cmd.plt, 1, 1) OF 
                           't': boite_pltz = box_h
                           ELSE : boite_pltz = [box_h, vertical_domain]
                        ENDCASE 
                     ENDIF 
                  ENDIF ELSE BEGIN
                     ; mask MSF as field=valmask
                     idx=where (fld eq valmask)
                     fld (idx) = !values.f_nan
                     boite_pltz = [box_h, vertical_domain]
                  ENDELSE 
               ENDIF ELSE BEGIN
                  CASE strmid(cmd.plt, 0, 2) OF
                  'xy': BEGIN & char = strmid(cmd.plt,3, 1) & ideb = 4 & END 
                  ELSE: BEGIN & char = strmid(cmd.plt,dimplot+1, 1) & ideb = dimplot+2 & END 
                  ENDCASE 
                  
                  CASE char OF
                     '#': BEGIN 
                                ; bathymetry mask for data
                        bat_name = strmid(cmd.plt, ideb, strlen(cmd.plt)-ideb)
                        full_name = grep('ls -al '+hom_idl+'grids | grep orca.'+bat_name+' | awk ''NR == 1 {print $0}''', ' ', 8)
                        full_name = strmid(full_name, 5, 50)
                        legz = full_name
                        legbox = full_name
                        IF cmd.var NE 'vozonbsf' THEN BEGIN 
                           maskzm = maskread(full_name, 'orca', /D3)
                           tmask = tmask*maskzm
                           IF dimplot EQ 1 THEN boite_pltz = box_h
                           IF dimplot EQ 2 THEN BEGIN
                              CASE strmid(cmd.plt, 1, 1) OF 
                                 't': boite_pltz = box_h
                                 ELSE : boite_pltz = [box_h, vertical_domain]
                              ENDCASE 
                           ENDIF 
                        ENDIF ELSE BEGIN
                                ; mask MSF as field=valmask
                           idx=where (fld eq valmask)
                           fld (idx) = !values.f_nan
                           boite_pltz = [box_h, vertical_domain]
                        ENDELSE  
                     END  
;                        '@': BEGIN
;                           ; density level
;                        END 
                     ELSE: BEGIN 
                                ; mask in box
                        zbox = def_box(cmd.plt, dimplot, legz, time_stride)
                        boite_pltz = zbox
                        IF dimplot EQ 2 THEN BEGIN
                           CASE strmid(cmd.plt, 1, 1) OF 
                              't': boite_pltz = zbox
                              ELSE : IF n_elements(boite_pltz) EQ 4 THEN boite_pltz = [zbox, vertical_domain]
                           ENDCASE 
                        ENDIF 
                        IF dimplot EQ 1 THEN BEGIN
                           CASE strmid(cmd.plt, 1, 1) OF 
                              't': boite_pltz = zbox
                              ELSE : BEGIN
                                 CASE strmid(cmd.plt, 0, 1) OF 
                                    'x': boite_pltz = zbox
                                    'y': boite_pltz = zbox
                                    ELSE: IF n_elements(boite_pltz) EQ 4 THEN boite_pltz = [zbox, vertical_domain]
                                 ENDCASE
                              END 
                           ENDCASE 
                        ENDIF 
                     END 
                  ENDCASE 
                  
              ENDELSE
               IF cmd.var EQ 'vozonbsf' AND msf_mean EQ 0 THEN BEGIN 
                   print, '     Mask_z special case for vozonbsf'
                   idx = where(gphit EQ max(gphit))
                   idx = idx MOD jpi
                   print, idx[0]
                   boite_pltz = [idx[0], idx[0], boite_pltz[2:5]]
                   domdef, boite_pltz, /xindex ; indice pour x
                   IF (size(fld))(0) EQ 3 THEN BEGIN
                       tmask[idx[0], firstyt:lastyt, *] = fld(*,*,1) LT valmask/10.
                   ENDIF ELSE BEGIN 
                       tmask[idx[0], firstyt:lastyt, *] = fld LT valmask/10.
                   ENDELSE  
               ENDIF
            END 
         ENDCASE 
      END   
      ELSE: BEGIN ; atmosphere or data
            
;         IF masked_data EQ 0 THEN tmask[*] = 1
         CASE strmid(cmd.plt, 0, 2) OF
            'xy': BEGIN & char = strmid(cmd.plt,3, 1) & ideb = 4 & END 
            ELSE: BEGIN & char = strmid(cmd.plt,dimplot+1, 1) & ideb = dimplot+2 & END 
         ENDCASE 


         IF char EQ '#' THEN BEGIN 
                                ; mask from file
            bat_name = strmid(cmd.plt, ideb, strlen(cmd.plt)-ideb)
            full_name = grep('ls -al '+hom_idl+'grids | grep '+cmd.grid+'.'+bat_name+' | awk ''NR == 1 {print $0}''', ' ', 8)
            full_name = strmid(full_name, strlen(cmd.grid)+1, 50)
            legz = full_name
            legbox = full_name
            maskzm = maskread(full_name, cmd.grid)
            maskzm = reverse(maskzm, 2)
            maskzm = reform(maskzm, jpi*jpj)
            tmask = reform(maskzm#replicate(1, jpk), jpi, jpj, jpk)
            boite_pltz = box_h
            IF dimplot EQ 2 THEN BEGIN
               CASE strmid(cmd.plt, 1, 1) OF 
                  't': boite_pltz = box_h
                  ELSE : boite_pltz = [box_h, hpa_min, hpa_max ]
               ENDCASE 
            ENDIF 
         ENDIF ELSE BEGIN 
           ; mask in box
            IF masked_data EQ 0 THEN BEGIN  
               CASE atmos_msk OF
               0: BEGIN & mskt = 'ALL points' & tmask[*] = 1 & END
               1: BEGIN & mskt = 'ocean points only' & END
               2: BEGIN & mskt = 'land points only' & tmask = 1-tmask & END
               ENDCASE 
               print, '    Using ',mskt, ' (atmos_msk =', atmos_msk, ')'
               print, ' '
            ENDIF  

            zbox = def_box(cmd.plt, dimplot, legz, time_stride)
            boite_pltz = zbox
            prof1 = gdept(0)
            prof2 = gdept(0)
            IF dimplot EQ 2 THEN BEGIN
               CASE strmid(cmd.plt, 1, 1) OF 
                  't': boite_pltz = zbox
                  ELSE : IF n_elements(boite_pltz) EQ 4 THEN boite_pltz = [zbox, hpa_min, hpa_max]
               ENDCASE 
            ENDIF 
         ENDELSE 
      END 
      ELSE:
   ENDCASE 
END