code: def isosurface(image, smoothing=0, threshold=None, connectivity=False):
      
      if smoothing:
          smImg = vtk.vtkImageGaussianSmooth()
          smImg.SetDimensionality(3)
          smImg.SetInputData(image)
          smImg.SetStandardDeviations(smoothing, smoothing, smoothing)
          smImg.Update()
          image = smImg.GetOutput()
  
      scrange = image.GetScalarRange()
      if scrange[1] > 1e10:
          print("Warning, high scalar range detected:", scrange)
  
      cf = vtk.vtkContourFilter()
      cf.SetInputData(image)
      cf.UseScalarTreeOn()
      cf.ComputeScalarsOn()
  
      if utils.isSequence(threshold):
          cf.SetNumberOfContours(len(threshold))
          for i, t in enumerate(threshold):
              cf.SetValue(i, t)
          cf.Update()
      else:
          if not threshold:
              threshold = (2 * scrange[0] + scrange[1]) / 3.0
          cf.SetValue(0, threshold)
          cf.Update()
  
      clp = vtk.vtkCleanPolyData()
      clp.SetInputConnection(cf.GetOutputPort())
      clp.Update()
      poly = clp.GetOutput()
  
      if connectivity:
          conn = vtk.vtkPolyDataConnectivityFilter()
          conn.SetExtractionModeToLargestRegion()
          conn.SetInputData(poly)
          conn.Update()
          poly = conn.GetOutput()
  
      a = Actor(poly, c=None).phong()
      a.mapper.SetScalarRange(scrange[0], scrange[1])
      return a
docs: Return a ``vtkActor`` isosurface extracted from a ``vtkImageData`` object.
  
      :param float smoothing: gaussian filter to smooth vtkImageData, in units of sigmas
      :param threshold:    value or list of values to draw the isosurface(s)
      :type threshold: float, list
      :param bool connectivity: if True only keeps the largest portion of the polydata
  
      .. hint:: |isosurfaces| |isosurfaces.py|_