How to find and split two surfaces made of roots of a determinant?

I have a matrix from which I want to find the roots for different values of [Delta]g0, [Delta]d0 and [Epsilon]. I know there will be a maximum of two roots for a set for a set of parameters and after that I would like to plot independently the surfaces made of the roots.

I have been using CountourPlot3D after some simplification on the determinant. The simplified determinant is Fusi.

In[202]:= FuSi = 
 FullSimplify[
  Re[Chop[FullSimplify[
     Chop [Det[
       Mraw[[Epsilon], [Delta]g0, [Delta]d0, U0, L0, 
        k0, [Kappa]00]]], {Element[[Delta]g0, Reals], 
      Element[[Delta]d0, Reals], 
      Element[[Epsilon], Reals], [Epsilon] > 0, [Epsilon] < 
       1, [Delta]g0 >= 0, [Delta]g0 <= 2*Pi, [Delta]d0 >= 
       0, [Delta]d0 <= 2*Pi}]]]]

Out[202]= Re[
 1/[Epsilon]^4 (5.0625 + 1. E^(-I ([Delta]d0 - [Delta]g0)) + 
    1. E^(I ([Delta]d0 - [Delta]g0)) + 
    1. E^(-I ([Delta]d0 + [Delta]g0)) + 
    1. E^(I ([Delta]d0 + [Delta]g0)) + 
    1. E^(I ([Delta]d0 - [Delta]g0) - 2. Sqrt[1 - [Epsilon]^2]) + 
    1. E^(-I [Delta]d0 + I [Delta]g0 - 2. Sqrt[1 - [Epsilon]^2]) + 
    1. E^(-I ([Delta]d0 + [Delta]g0) - 2. Sqrt[1 - [Epsilon]^2]) + 
    1. E^(I ([Delta]d0 + [Delta]g0) - 2. Sqrt[1 - [Epsilon]^2]) - 
    2. E^(-I ([Delta]d0 + [Delta]g0) - 1. Sqrt[1 - [Epsilon]^2]) - 
    2. E^(I ([Delta]d0 + [Delta]g0) - 1. Sqrt[1 - [Epsilon]^2]) - 
    19.125 [Epsilon]^2 + 
    4. E^(-I [Delta]d0 - 1. Sqrt[1 - [Epsilon]^2]) [Epsilon]^2 + 
    4. E^(I [Delta]d0 - 1. Sqrt[1 - [Epsilon]^2]) [Epsilon]^2 + 
    4. E^(-I [Delta]g0 - 1. Sqrt[1 - [Epsilon]^2]) [Epsilon]^2 + 
    4. E^(I [Delta]g0 - 1. Sqrt[1 - [Epsilon]^2]) [Epsilon]^2 + 
    18.0625 [Epsilon]^4 + 
    E^(-I [Delta]g0) (2.25 - 4.25 [Epsilon]^2) + 
    E^(I [Delta]g0) (2.25 - 4.25 [Epsilon]^2) + 
    E^(I [Delta]d0) (2.25 - 4.25 [Epsilon]^2) + 
    E^(-I [Delta]d0) (2.25 - 4.25 [Epsilon]^2) + 
    E^(-I [Delta]d0 + I [Delta]g0 - 
      1. Sqrt[1 - [Epsilon]^2]) (2. - 4. [Epsilon]^2) + 
    E^(I ([Delta]d0 - [Delta]g0) - 
      1. Sqrt[1 - [Epsilon]^2]) (2. - 4. [Epsilon]^2) + 
    E^(I [Delta]d0 - 
      2. Sqrt[1 - [Epsilon]^2]) (-2.25 + 0.25 [Epsilon]^2) + 
    E^(-I [Delta]d0 - 
      2. Sqrt[1 - [Epsilon]^2]) (-2.25 + 0.25 [Epsilon]^2) + 
    E^(-I [Delta]g0 - 
      2. Sqrt[1 - [Epsilon]^2]) (-2.25 + 0.25 [Epsilon]^2) + 
    E^(I [Delta]g0 - 
      2. Sqrt[1 - [Epsilon]^2]) (-2.25 + 0.25 [Epsilon]^2) + 
    E^(-2. Sqrt[
      1 - [Epsilon]^2]) (5.0625 - 1.125 [Epsilon]^2 + 
       0.0625 [Epsilon]^4) + 
    E^(-1. Sqrt[
      1 - [Epsilon]^2]) (1.875 - 11.75 [Epsilon]^2 + 
       1.875 [Epsilon]^4))]

From this simplication I can extract the roots using ContourPlot3D on a set of parameters.

In[215]:= Det6[phil_?NumericQ, phir_?NumericQ, energy_?NumericQ] := 
 FuSi /. {[Delta]g0 -> phil, [Delta]d0 -> phir, [Epsilon] -> 
    energy}

In[216]:=spect= ContourPlot3D[Det6[phil,phir,energy]==0,{phil,0,2*Pi},{phir,0,2*Pi},{energy,0.01,0.99},Contours->0]

Which gives :

So far I have the surface define by the roots and now I want to extract a the set of coordinates corresponding to each of the two surfaces independently. I have tried using

data = First@
  Cases[spect, GraphicsComplex[points_, ___] :> points, Infinity]  

But it gives a mixed list of points from which I can’t split the two surfaces when I try to filter it using:

In[440]:=sorteddata = SortBy[data, {#[[1]] &, #[[2]] &, #[[3]] &}]
In[441]:=BotSurface= 
 DeleteDuplicates[
  sorteddata, (#1[[1]] == #2[[1]] && #1[[2]] == #2[[2]] && #1[[3]] < 
#2[[3]]) &]

The problem is that my dataset doesn’t have the same XY points for each point on the surfaces. Some points with specific XY coordinates calculated on the top surface are not calculated on the bottom surface. So I can’t simply discriminate points like that.

I’m pretty sure there is a function in Mathematica to help me do that but I can’t figure out what is it ?
I’m not against a completely other method also.

EDIT:
Adding some definitions if you want to run the first In[202]

In[200]:= u[[Epsilon]_, [Delta]_] = (E^(I [Delta]/2) Sqrt[
      1 + Sqrt[1 - 1/[Epsilon]^2]])/Sqrt[2];
v[[Epsilon]_, [Delta]_] = (E^(-I [Delta]/2) Sqrt[
      1 - Sqrt[1 - 1/[Epsilon]^2]])/Sqrt[2];
[Kappa][[Epsilon]_, [Kappa]0_] = [Kappa]0*Sqrt[[Epsilon]^2 - 1];
[Delta]d0 =.;
[Delta]g0 =.;
U0 = 0.25;
L0 = 0.5;
[Kappa]00 = 1;
[Epsilon]0 =.;
x0 = -L0/2;
k0 = (2 [Pi]/10^-4 + [Pi]/2*0)/L0;
phil =.;
phir =.;
energy =.;

In[201]:= Mraw[[Epsilon]_, [Delta]g_, [Delta]d_, U_, L_, k_, [Kappa]0_] = ({
    {u[[Epsilon], [Delta]g]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], 
     v[[Epsilon], -[Delta]g]*
      Exp[-I*L/
         2*(k - [Kappa][[Epsilon], [Kappa]0])], -u[[Epsilon], 0]*
      Exp[-I*L/
         2*(k + [Kappa][[Epsilon], [Kappa]0])], -u[[Epsilon], 0]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], -v[[Epsilon],
         0]*Exp[-I*
        L/2*(k - [Kappa][[Epsilon], [Kappa]0])], -v[[Epsilon], 0]*
      Exp[I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], 0, 0},
    {v[[Epsilon], [Delta]g]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], 
     u[[Epsilon], -[Delta]g]*
      Exp[-I*L/
         2*(k - [Kappa][[Epsilon], [Kappa]0])], -v[[Epsilon], 0]*
      Exp[-I*L/
         2*(k + [Kappa][[Epsilon], [Kappa]0])], -v[[Epsilon], 0]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], -u[[Epsilon],
         0]*Exp[-I*
        L/2*(k - [Kappa][[Epsilon], [Kappa]0])], -u[[Epsilon], 0]*
      Exp[I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], 0, 0},
    {0, 0, -u[[Epsilon], 0]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], -u[[Epsilon],
         0]*Exp[-I*
        L/2*(k + [Kappa][[Epsilon], [Kappa]0])], -v[[Epsilon], 0]*
      Exp[I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], -v[[Epsilon],
         0]*Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], 
     u[[Epsilon], [Delta]d]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], 
     v[[Epsilon], -[Delta]d]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])]},
    {0, 0, -v[[Epsilon], 0]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], -v[[Epsilon],
         0]*Exp[-I*
        L/2*(k + [Kappa][[Epsilon], [Kappa]0])], -u[[Epsilon], 0]*
      Exp[I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], -u[[Epsilon],
         0]*Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], 
     v[[Epsilon], [Delta]d]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], 
     u[[Epsilon], -[Delta]d]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])]},
    {(1 + I*U)*u[[Epsilon], [Delta]g]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      v[[Epsilon], -[Delta]g]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      u[[Epsilon], 0]*
      Exp[-I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      u[[Epsilon], 0]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      v[[Epsilon], 0]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      v[[Epsilon], 0]*
      Exp[I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], 0, 0},
    {(1 + I*U)*v[[Epsilon], [Delta]g]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      u[[Epsilon], -[Delta]g]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      v[[Epsilon], 0]*
      Exp[-I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      v[[Epsilon], 0]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      u[[Epsilon], 0]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      u[[Epsilon], 0]*
      Exp[I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], 0, 0},
    {0, 0, (-1 + I*U)*u[[Epsilon], 0]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      u[[Epsilon], 0]*
      Exp[-I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      v[[Epsilon], 0]*
      Exp[I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      v[[Epsilon], 0]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      u[[Epsilon], [Delta]d]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      v[[Epsilon], -[Delta]d]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])]},
    {0, 0, (-1 + I*U)*v[[Epsilon], 0]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      v[[Epsilon], 0]*
      Exp[-I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      u[[Epsilon], 0]*
      Exp[I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      u[[Epsilon], 0]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])], (1 + I*U)*
      v[[Epsilon], [Delta]d]*
      Exp[I*L/2*(k + [Kappa][[Epsilon], [Kappa]0])], (-1 + I*U)*
      u[[Epsilon], -[Delta]d]*
      Exp[-I*L/2*(k - [Kappa][[Epsilon], [Kappa]0])]}
   });