code: def pip(self, points, sorted_col=0, radius=0):
          
          xy = points[:, :2]
          n_points = xy.shape[0]
          index = np.arange(n_points, dtype = int)
          b = self.domain
          b[0] = b[0] - radius
          b[1] = b[1] + radius
          
          # Slicing the sorted column
          k = np.searchsorted(xy[:, sorted_col],
                              (b[0, sorted_col], b[1, sorted_col]+1e-10))
          xy = xy[k[0]:k[1]]
          index = index[k[0]:k[1]]
          
          # solution
          k = index[self.path.contains_points(xy, radius=radius)]
          sol = np.zeros(n_points, dtype=bool)
          sol[k] = True
                  
          return sol
docs: Point-in-Polygon for the z=0 projection. This function enhances
          the performance of ``Polygon.contains()`` by verifying only the 
          points which are inside the bounding box of the polygon. To do 
          it fast, it needs the points array to be already sorted by one
          column.
          
          :param points: list of *(x, y, z) or (x, y)* coordinates of the
              points to check. (The z value will not be taken into 
              account).
          :type points: ndarray (shape=(N, 2 or 3))
          :param sorted_col: Index of the sorted column (0 or 1).
          :type sorted_col: int
          :param radius: Enlarge Polygons domain by a specified quantity.
          :type radius: float
          :returns: Which points are inside the polygon.
          :rtype: ndarray (dtpye=bool)
          
          .. warning:: By default pip considers that the set of points is
              currently sorted by the first column.
          .. warning:: This method only works if the polygon has been 
              locked (:func:`lock`).