;------------------------------------------------------------
;------------------------------------------------------------
;------------------------------------------------------------
;+
; NAME: maskread
;   
;
; PURPOSE:
;       lit un fichier bathymetry ORCA ou un mask atmos
;
; CATEGORY:
;       lecture de fichier
;
; INPUTS:

; KEYWORD PARAMETERS:

; COMMON BLOCKS:
;       common.pro
;
; SIDE EFFECTS:Retourne -1 en cas d'erreur.
;
; RESTRICTIONS:
;
; EXAMPLE:
;
; MODIFICATION HISTORY: Maurice Imbard (mimlod@ipsl.jussieu.fr)
;			12/05/99
;-
;--------------------------------------------------------------
;--------------------------------------------------------------
;--------------------------------------------------------------
function maskread, bathy, model, D3 = d3
;;
@common
@com_eg
;; 
;;
   IF n_elements(nbathys) EQ 0 THEN BEGIN
      read_bathy = 'yes'
   ENDIF ELSE BEGIN 
      i = 0
      WHILE i LE nbathys-1 DO BEGIN
         IF bathy_read[i].name EQ bathy THEN BEGIN
            zbat = bathy_read[i].mask
            IF keyword_set(D3) THEN BEGIN
               zbat = 0 > zbat 
               ; Build 3D from pseudo 3D
               zbat=reform(zbat[*]#(1./(indgen(jpkglo)+1)),(size(zbat))[1],(size(zbat))[2],jpkglo)
               zbat=floor(zbat)<1
            ENDIF ELSE BEGIN 
               zbat = zbat <  1
            ENDELSE 
            return, zbat
         ENDIF ELSE BEGIN 
            read_bathy = 'yes'
         ENDELSE
         i = i+1
      ENDWHILE 
   ENDELSE 

   IF read_bathy EQ 'yes' THEN BEGIN 
      z = -1
;;
;---------------------------------------------------------------------
; definition de la taille du fichier a aller chercher: jpidta,jpjdta,jpkdta...
;---------------------------------------------------------------------
      if n_elements(jpidta) EQ 0 THEN BEGIN 
         if n_elements(ixmindta) EQ 0 OR n_elements(ixmaxdta) EQ 0 then $
          jpidta = jpiglo else jpidta = ixmaxdta-ixmindta+1
      endif
      if n_elements(jpjdta) EQ 0 THEN BEGIN 
         if n_elements(iymindta) EQ 0 OR n_elements(iymaxdta) EQ 0 then $
          jpjdta = jpjglo else jpjdta = iymaxdta-iymindta+1
      endif
      if n_elements(jpkdta) EQ 0 THEN BEGIN 
         if n_elements(izmindta) EQ 0 OR n_elements(izmaxdta) EQ 0 then $
          jpkdta = jpkglo else jpkdta = izmaxdta-izmindta+1
      endif
;--------------------------------------------------------------------
; ouverture du fichier a l'adresse s_fichier
;--------------------------------------------------------------------
;;
      iname = bathy
      iname_file = hom_idl+'grids/'+model+'.'+iname
      openr, iunit,  /get_lun,iname_file
      print, '    Reading ', model+'.'+iname, ' mask '
      print, ' '
;---------------------------------------------------------------------
; lecture des champs directement vers le champ
;---------------------------------------------------------------------
      ifreq = 40
      ifin  = jpiglo/ifreq+1
      irest = jpiglo-(ifin-1)*ifreq
      zbati  = intarr(ifreq)
      zbati2 = intarr(irest)
      ;;
      zbat    = intarr(jpiglo,jpjglo)
;
      readf, iunit, FORMAT = '(16x,2i8)', iim, ijm
      readf, iunit, FORMAT = '(/)'
      readf, iunit, FORMAT = '(/)'
      il1 = 0
      FOR jn = 1, ifin-1 DO BEGIN
         readf, iunit, FORMAT = '(/)'
         readf, iunit, FORMAT = '(/)'
         il2 = min([ jpiglo-1, il1+ifreq-1] )
         readf, iunit, FORMAT = '(/)'
         readf, iunit, FORMAT = '(/)'
         readf, iunit, FORMAT = '(/)'
         il3 = il2-(jn-1)*ifreq
         iformat = string('(', il3+2, 'i3)')
         FOR jj =  jpjglo-1, 0, -1  DO BEGIN
            zbati(0:ifreq-1) = 0
            readf, iunit, FORMAT = iformat, ij, zbati
            zbat(il1:il2, jj) = zbati(0:il3)
         END
         il1 = il1 + ifreq
      END
      readf, iunit, FORMAT = '(/)'
      readf, iunit, FORMAT = '(/)'
      il2 = min([ jpiglo-1, il1+ifreq-1] )
      readf, iunit, FORMAT = '(/)'
      readf, iunit, FORMAT = '(/)'
      readf, iunit, FORMAT = '(/)'
      il3 = il2-(ifin-1)*ifreq
      iformat = string('(', il3+2, 'i3)')
      FOR jj =  jpjglo-1, 0, -1  DO BEGIN
         zbati2(0:irest-1) = 0
         readf, iunit, FORMAT = iformat, ij, zbati2
         zbat(il1:il2, jj) = zbati2(0:il3)
      END
      ;;
;---------------------------------------------------------------------
; on initialise les ixmindta, iymindta  au besoin
;---------------------------------------------------------------------
;;
      ixmindta = 0
      ixmaxdta = jpiglo-1
      iymindta = 0
      iymaxdta = jpjglo-1
      izmindta = 0
      izmaxdta = jpkglo-1
;---------------------------------------------------------------------
; on reduit z selon les valeurs de ixmindta, iymindta, ...+ shift 
;---------------------------------------------------------------------
      z = zbat[ixminmesh-ixmindta:ixmaxmesh-ixmindta $
               ,iyminmesh-iymindta:iymaxmesh-iymindta ]
      z = shift(z,key_shift, 0)
;---------------------------------------------------------------------
      IF n_elements(nbathys) EQ 0 THEN BEGIN
         nbathys = 1
         a = z*0
         bathy_read = {name: '', mask:a}
         bathy_read = replicate(bathy_read, 4)
         bathy_read[0].name = bathy
         bathy_read[0].mask = z
      ENDIF ELSE BEGIN 
         IF nbathys LE 3 THEN BEGIN
            nbathys = nbathys + 1
         ENDIF ELSE BEGIN 
            nbathys = 1
         ENDELSE 
         bathy_read[nbathys-1].name = bathy
         bathy_read[nbathys-1].mask = z 
      ENDELSE 

;---------------------------------------------------------------------
fini:
      free_lun,iunit
      close, iunit
;---------------------------------------------------------------------
sortie:
      zbat = z
      IF keyword_set(D3) THEN BEGIN
         zbat = 0 > zbat 
         ; Build 3D from pseudo 3D
         zbat=reform(zbat[*]#(1./(indgen(jpkglo)+1)),(size(zbat))[1],(size(zbat))[2], jpkglo)
         zbat=floor(zbat)<1
      ENDIF ELSE BEGIN 
         zbat = zbat <  1
      ENDELSE 
      return, zbat
   ENDIF 
END