;+
    ;
    ;   GEOST_1D : Compute geostrophic speeds from a sea level dataset <br />
    ;
    ;   Reference : Powell, B. S., et R. R. Leben (2004), An Optimal Filter for <br />
    ;   Geostrophic Mesoscale Currents from Along-Track Satellite Altimetry, <br />
    ;   Journal of Atmospheric and Oceanic Technology, 21(10), 1633-1642.
    ;
    ; @param lon {in}{optional}{type:NUMERIC} longitude in degrees
    ; @param lat {in}{optional}{type:NUMERIC} latitude in degrees
    ; @param dst {in}{optional}{type:NUMERIC} along-track distance.
    ; @param z {in}{required}{type:NUMERIC} sea level surface. Can either be absolute<br />
    ;          values (SSH) or relative values (SLA). This MUST be given in METERS.
    ; @keyword strict {in}{optional}{type:BOOLEAN} If True, compute gradient at mid-distance.
    ; @keyword pl04 {in}{optional}{type:BOOLEAN} If True, use the Powell & Leben 2004 method.
    ;          
    ; @returns Geostrophic velocity component, positive eastward
    ;
    ;
    ;
    ; @author Renaud DUSSURGET, LEGOS/CTOH
    ; @history Created Sep. 2009 from genweights.m (Brian Powell (c) 2004, <br />
    ;   University of Colorado, Boulder)<br />
    ;   Modified May 2010 to be compliant with 20Hz datasets (p & n can vary).<br />
    ;     Warining may also be issued for data with holes within the width of the<br />
    ;     window.<br />
    ;   Modified June 2010 to include filtering window width in KM instead of nb. of<br />
    ;     points (Equivalent btw. 1Hz and 20Hz data).<br />
    ;
    ; @uses CALCUL_DISTANCE, EXIST, GENWEIGTHS, SETINTERSECTION, SETUNION, <br />
    ;   OPTIMAL_SLOPE, GRAVITY, CORIOLIS, TRACK_ORIENT
    ;
    ; @example dummy1=geost_1D(lon,lat,sla,pl04=True,p=11,q=11) :<br />
    ;   Return along-track velocity anomalies using a 11km by 11km Powell & Leben 2004 filter window <br />
    ;          dummy2=geost_1D(dst,sla,strict=True) :<br />
    ;   Return along-track velocity anomalies computed at mid-distance <br />
    ;
    ;-

def geost_1d(*args,**kwargs) : #(lon,lat,nu): OR (dst,nu)
    """
    ;+
    ;
    ;   GEOST_1D : Compute geostrophic speeds from a sea level dataset <br />
    ;
    ;   Reference : Powell, B. S., et R. R. Leben (2004), An Optimal Filter for <br />
    ;   Geostrophic Mesoscale Currents from Along-Track Satellite Altimetry, <br />
    ;   Journal of Atmospheric and Oceanic Technology, 21(10), 1633-1642.
    ;
    ; @param lon {in}{optional}{type:NUMERIC} longitude in degrees
    ; @param lat {in}{optional}{type:NUMERIC} latitude in degrees
    ; @param dst {in}{optional}{type:NUMERIC} along-track distance.
    ; @param z {in}{required}{type:NUMERIC} sea level surface. Can either be absolute<br />
    ;          values (SSH) or relative values (SLA). This MUST be given in METERS.
    ; @keyword strict {in}{optional}{type:BOOLEAN} If True, compute gradient at mid-distance.
    ; @keyword pl04 {in}{optional}{type:BOOLEAN} If True, use the Powell & Leben 2004 method.
    ;          
    ; @returns Geostrophic velocity component, positive eastward
    ;
    ;
    ;
    ; @author Renaud DUSSURGET, LEGOS/CTOH
    ; @history Created Sep. 2009 from genweights.m (Brian Powell (c) 2004, <br />
    ;   University of Colorado, Boulder)<br />
    ;   Modified May 2010 to be compliant with 20Hz datasets (p & n can vary).<br />
    ;     Warining may also be issued for data with holes within the width of the<br />
    ;     window.<br />
    ;   Modified June 2010 to include filtering window width in KM instead of nb. of<br />
    ;     points (Equivalent btw. 1Hz and 20Hz data).<br />
    ;
    ; @uses CALCUL_DISTANCE, EXIST, GENWEIGTHS, SETINTERSECTION, SETUNION, <br />
    ;   OPTIMAL_SLOPE, GRAVITY, CORIOLIS, TRACK_ORIENT
    ;
    ; @example dummy1=geost_1D(lon,lat,sla,pl04=True,p=11,q=11) :<br />
    ;   Return along-track velocity anomalies using a 11km by 11km Powell & Leben 2004 filter window <br />
    ;          dummy2=geost_1D(dst,sla,strict=True) :<br />
    ;   Return along-track velocity anomalies computed at mid-distance <br />
    ;
    ;-
    """
    
    lon = args[0]
    lat = args[1]
    
    dst = args[2] if len(args) == 4 else calcul_distance(lat,lon) * 1e3 #distance in meters
    nu = args [3] if len(args) == 4 else args[2]
    
    isVector = len(np.shape(nu)) == 1 
    
    #Reshape nu if vector
    if isVector : nu=np.reshape(nu,(len(nu),1))
    
    nt = np.shape(nu)[1] if not isVector else 1
    sh = nu.shape
    nufilt=np.ma.array(np.empty(sh),mask=True,dtype=nu.dtype)
    
    pl04 = kwargs.pop('pl04',False)
    filter = kwargs.pop('filter', None)
    strict = kwargs.pop('strict',False)
    verbose = kwargs.pop('verbose',False)
    
    if filter is not None :
        for t in np.arange(nt) : 
            nufilt[:,t] =loess(nu[:,t],dst,filter*1e3)
        nu=nufilt
    
    if pl04 :
        ug = np.ma.array(np.empty(sh),mask=True,dtype=nu.dtype)
        for t in np.arange(nt) :
            ug[:,t] = powell_leben_filter_km(lon,lat,nu[:,t],verbose=verbose,**kwargs)
        if isVector : ug=ug.flatten()
        return ug
    
    #If strict option is set to True, compute gradients at mid-distance between points
    if strict :
        lon = (lon[1:] - lon[:-1])/2. + lon[0:-1]
        lat = (lat[1:] - lat[:-1])/2. + lat[0:-1]
    
    
    #Compute gravitational & coriolis forces
    if strict : sh = (sh[0]-1,sh[1])
    g = np.repeat(gravity(lat),nt).reshape(sh)
    f = np.repeat(coriolis(lat),nt).reshape(sh)
    
    #Compute SSH 1st derivative
#    dh = deriv(dst,nu) #(deriv is very bad...)
    
    
    dh = np.ma.array(np.empty(sh),mask=True,dtype=nu.dtype)
    for t in np.arange(nt) :
        dh[:,t] = (nu[1:,t] - nu[:-1,t])/(dst[1:] - dst[:-1]) if strict else deriv(dst,nu[:,t])
      
    #Compute geostrophy
#    print f
#    print g
#    print dh
    ug = - (g*dh) / (f)
  
    #Inverse sign of ug for descending tracks as Coriolis is oriented to the right
    #northward
    if (not track_orient(lon,lat)) : #descending tracks
        ug *=-1
    
    if isVector : ug=ug.flatten()
    
    return (lon,lat,ug) if strict else ug