;+
;
; @file_comments
;
; @categories
; Grid
;
; @param XAXIS
;
; @param YAXIS
;
; @param ZAXIS
;
; @param LINETYPE
;
; @keyword FORTHEMASK
;
; @keyword WPOINT
;
; @keyword FIRSTS
;
; @keyword LASTS
;
; @keyword PTTYPE
;
; @uses
; common.pro
;
; @history
; Sebastien Masson (smasson\@lodyc.jussieu.fr)
;                      2001-06
;
; @version
; $Id$
;
; @todo seb
;
;-
;
PRO computehopegrid, xaxis, yaxis, zaxis, linetype, FORTHEMASK = forthemask, WPOINT = wpoint, FIRSTS = firsts, LASTS = lasts, PTTYPE = pttype
;---------------------------------------------------------
;
  compile_opt idl2, strictarrsubs
;
@cm_4mesh
@cm_4data
  IF NOT keyword_set(key_forgetold) THEN BEGIN
@updatenew
  ENDIF
;---------------------------------------------------------
;------------------------------------------------------------
   if not (keyword_set(scalar)*keyword_set(vector)) then scalar=1
;------------------------------------------------------------
   jpiglo = n_elements(xaxis)
   jpjglo = n_elements(yaxis)
   jpkglo = n_elements(zaxis)
;
   jpi = jpiglo
   jpj = jpjglo
   jpk = jpkglo
;
   if NOT keyword_set(firsts) then firsts = [0, 0, 0]
   if NOT keyword_set(lasts) then lasts = [jpi-1, jpj-1, jpk-1]
;
; determination of the grid type and of the point type
;
   if keyword_set(pttype) then vargrid = pttype
   if linetype EQ 'odd-even' then key_gridtype = 'e' ELSE key_gridtype = 'c'
;
; computation of the horizontal grid
;
   if key_gridtype EQ 'e' then begin
      if vargrid EQ 'T' then begin
         glamt=xaxis
         glamt = glamt#replicate(1, jpj)
         xstep=(glamt[1,0]-glamt[0,0])/2.
         glamt[*,2*lindgen((jpj)/2)]=glamt[*,2*lindgen(jpj/2)]-xstep
         glamu=glamt+xstep
      ENDIF ELSE BEGIN
         glamu=xaxis
         glamu = glamu#replicate(1, jpj)
         xstep=(glamu[1,0]-glamu[0,0])/2.
         glamu[*,2*lindgen((jpj)/2)]=glamu[*,2*lindgen(jpj/2)]-xstep
         glamt=glamu-xstep
      ENDELSE
;zoom
      glamt = glamt[firsts[0]:lasts[0], firsts[1]:lasts[1]]
      glamu = glamu[firsts[0]:lasts[0], firsts[1]:lasts[1]]
      jpiglo = lasts[0]-firsts[0]+1
      jpi = jpiglo
      jpjglo = lasts[1]-firsts[1]+1
      jpj = jpjglo
      glamv = glamu
      glamf = glamu
      gphit = yaxis[firsts[1]:lasts[1]]
      gphit = replicate(1, jpi)#gphit
      gphif = gphit
      gphiu = gphit
      gphiv = gphif
   ENDIF ELSE BEGIN
      if vargrid eq 'T' then begin
         glamt=xaxis
         glamt = glamt#replicate(1, jpj)
         glamu=glamt+(glamt[1,0]-glamt[0,0])/2.
      ENDIF ELSE BEGIN
         glamu=xaxis
         glamu = glamu#replicate(1, jpj)
         xstep=(glamu[1,0]-glamu[0,0])/2.
         glamt=glamu-(glamu[1,0]-glamu[0,0])/2.
      ENDELSE
;zoom
      glamt = glamt[firsts[0]:lasts[0], firsts[1]:lasts[1]]
      glamu = glamu[firsts[0]:lasts[0], firsts[1]:lasts[1]]
      jpiglo = lasts[0]-firsts[0]+1
      jpi = jpiglo
      jpjglo = lasts[1]-firsts[1]+1
      jpj = jpjglo
      glamv = glamt
      glamf = glamu
      gphit = yaxis[firsts[1]:lasts[1]]
      gphit = replicate(1, jpi)#gphit
      gphiu = gphit
      if jpj GT 1 then begin
         gphiv = gphit+(gphit[0, 1]-gphit[0, 0])/2.
         gphif = gphit+(gphit[0, 1]-gphit[0, 0])/2.
      ENDIF ELSE BEGIN
         gphiv = gphit
         gphif = gphit
      ENDELSE
   ENDELSE
;
; computation of the vertical grid
;
   if keyword_set(wpoint) then begin
      gdepw = zaxis
      if jpk ne 1 then begin
         gdept=(shift(gdepw, -1)+gdepw)/2.
         gdept[jpk-1]=gdepw[jpk-1]+.5*(gdepw[jpk-1]-gdepw[jpk-2])
      endif else gdept=zaxis
   endif else begin
      gdept = zaxis
      if jpk ne 1 then begin
         gdepw=(shift(gdept, 1)+gdept)/2.
         gdepw[0]=0
      endif else gdepw = zaxis
   endelse
;
; computation of the vertical scale factors
;
   if jpk ne 1 then begin
      e3t = abs(shift(gdepw,-1)-gdepw)
      e3t[jpk-1] = abs(gdept[jpk-1]-gdepw[jpk-1])
      e3w = abs(gdept-shift(gdept,1))
      e3w[0] = abs(gdept[0]-gdepw[0])
   endif else begin
      e3t=1
      e3w=1
   endelse
; zoom
   gdept = gdept[firsts[2]:lasts[2]]
   gdepw = gdepw[firsts[2]:lasts[2]]
   e3t = e3t[firsts[2]:lasts[2]]
   e3w = e3w[firsts[2]:lasts[2]]
   jpkglo = lasts[2]-firsts[2]+1
   jpk = jpkglo
;
; computation of the horizontal scale factors
;
   e1t = replicate(1b,jpi,jpj)
   e2t = replicate(1b,jpi,jpj)
   e1u = e1t
   e2u = e2t
   e1v = e1t
   e2v = e2t
   e1f = e1t
   e2f = e2t
;
; mask
;
   tmask = replicate(1b, jpi, jpj, jpk)
   if keyword_set(forthemask) then BEGIN
      land=where(forthemask ge valmask/10)
      if land[0] ne -1 then tmask[land] = 0b
   endif
   umaskred = replicate(1, jpj, jpk)
   vmaskred = replicate(1, jpi, jpk)
   fmaskredy = replicate(1, jpj, jpk)
   fmaskredx = replicate(1, jpi, jpk)
   @updateold
   domdef
   if keyword_set(firsts) AND keyword_set(lasts) then BEGIN
      if vargrid EQ 'T' then BEGIN
         if jpj GT 1 then begin
            lon1 = min(glamt[0, [0, 1]])
            lon2 = max(glamt[jpi-1, [0, 1]])
         endif ELSE BEGIN
            lon1 = min(glamt[0, 0])
            lon2 = max(glamt[jpi-1, [0]])
         ENDELSE
      ENDIF ELSE BEGIN
         if jpj GT 1 then begin
            lon1 = min(glamu[0, [0, 1]])
            lon2 = max(glamu[jpi-1, [0, 1]])
         endif ELSE BEGIN
            lon1 = min(glamu[0, 0])
            lon2 = max(glamu[jpi-1, 0])
         ENDELSE
      ENDELSE
      lat1 = min([gphit[0, 0], gphit[0, jpj-1]])
      lat2 = max([gphit[0, 0], gphit[0, jpj-1]])
      domdef, lon1, lon2, lat1, lat2, gdepw[0], gdept[jpk-1], gridtype = vargrid
   ENDIF
;
   ixminmesh  =0l
   ixmaxmesh  =long(jpi-1)
   iyminmesh  =0l
   iymaxmesh  =long(jpj-1)
   izminmesh  =0l
   izmaxmesh  =long(jpk-1)
;----------------------------------------------------------
; for the triangulation
;----------------------------------------------------------
   key_periodic = glamt[0] EQ ((glamt[jpi-1]+glamt[+1]-glamt[0]) MOD 360)
   if jpi gt 4 AND jpj GT 4 then begin
      triangles_list = triangule(shifted = glamt[0, 0] LT glamt[0, 1])
      twin_corners_up=-1
      twin_corners_dn=-1
   ENDIF ELSE BEGIN
      triangles_list=-1
      twin_corners_up=-1
      twin_corners_dn=-1
   ENDELSE
;------------------------------------------------------------
  IF NOT keyword_set(key_forgetold) THEN BEGIN
   @updateold
  ENDIF
;------------------------------------------------------------
   return
end