docs: def calc_kInkOut_IsoFlux(self, lPoly, lVIn=None, Lim=None,
                               kInOut=True):
          """ Calculate the intersection points of each ray with each isoflux
  
          The isofluxes are provided as a list of 2D closed polygons
  
          The intersections are the inward and outward intersections
          They are retruned as two np.ndarrays: kIn and kOut
          Each array contains the length parameter along the ray for each isoflux
  
          Parameters
          ----------
  
  
          Returns
          -------
  
          """
  
          # Preformat input
          nPoly, lPoly, lVIn = self._kInOut_IsoFlux_inputs_usr(lPoly, lVIn=lVIn)
  
          # Prepare output
          kIn = np.full((self.nRays,nPoly), np.nan)
          kOut = np.full((self.nRays,nPoly), np.nan)
  
          # Compute intersections
          assert(self._method in ['ref', 'optimized'])
          if self._method=='ref':
              for ii in range(0,nPoly):
                  largs, dkwd = self._kInOut_IsoFlux_inputs([lPoly[ii]],
                                                            lVIn=[lVIn[ii]])
                  out = _GG.SLOW_LOS_Calc_PInOut_VesStruct(*largs, **dkwd)
                  PIn, POut, kin, kout, VperpIn, vperp, IIn, indout = out
                  kIn[:,ii], kOut[:,ii] = kin, kout
          elif self._method=="optimized":
              for ii in range(0,nPoly):
                  largs, dkwd = self._kInOut_IsoFlux_inputs([lPoly[ii]],
                                                            lVIn=[lVIn[ii]])
  
                  out = _GG.LOS_Calc_PInOut_VesStruct(*largs, **dkwd)
                  kin, kout, _, _ = out
                  kIn[:,ii], kOut[:,ii] = kin, kout
          if kInOut:
              indok = ~np.isnan(kIn)
              ind = np.zeros((self.nRays,nPoly), dtype=bool)
              kInref = np.tile(self.kIn[:,np.newaxis],nPoly)
              kOutref = np.tile(self.kOut[:,np.newaxis],nPoly)
              ind[indok] = (kIn[indok]<kInref[indok]) | (kIn[indok]>kOutref[indok])
              kIn[ind] = np.nan
  
              ind[:] = False
              indok[:] = ~np.isnan(kOut)
              ind[indok] = (kOut[indok]<kInref[indok]) | (kOut[indok]>kOutref[indok])
              kOut[ind] = np.nan
  
          return kIn, kOut