code: def gen_topomask(h, lon, lat, dx=1.0, kind="linear", plot=False):
      
  
      import gsw
      from scipy.interpolate import interp1d
  
      h, lon, lat = list(map(np.asanyarray, (h, lon, lat)))
      # Distance in km.
      x = np.append(0, np.cumsum(gsw.distance(lon, lat)[0] / 1e3))
      h = -gsw.z_from_p(h, lat.mean())
      Ih = interp1d(x, h, kind=kind, bounds_error=False, fill_value=h[-1])
      xm = np.arange(0, x.max() + dx, dx)
      hm = Ih(xm)
  
      return xm, hm
docs: Generates a topography mask from an oceanographic transect taking the
      deepest CTD scan as the depth of each station.
  
      Inputs
      ------
      h : array
          Pressure of the deepest CTD scan for each station [dbar].
      lons : array
             Longitude of each station [decimal degrees east].
      lat : Latitude of each station. [decimal degrees north].
      dx : float
           Horizontal resolution of the output arrays [km].
      kind : string, optional
             Type of the interpolation to be performed.
             See scipy.interpolate.interp1d documentation for details.
      plot : bool
             Whether to plot mask for visualization.
  
      Outputs
      -------
      xm : array
           Horizontal distances [km].
      hm : array
           Local depth [m].
  
      Author
      ------
      André Palóczy Filho (paloczy@gmail.com) --  October/2012