I am using the following code to find the eigenvalues of a tightbinding matrix

ClearAll["Global`*"] (*Clear all variables*)
h := 6.63*10^-34; (* J.s, Planck's constant *)
hbar := (663*10^-34)/(100*2 \[Pi]); (* J.s, Reduced Planck's constant *)
q := 16/10*10^-19 ;  (* Coulomb, electron charge *)
m0 := 91/10*10^-31; (* kg free electron rest mass *)
a := 5.43*10^-10;
vss\[Sigma] := 2.03;
vsp\[Sigma] := 2.55;
vpp\[Sigma] := 4.55;
vpp\[Pi] := 1.09;
esa := -13.55;
esb := -13.55;
epa := -6.52;
epb := -6.52;
n1 := a/4*{1, 1, 1};
n2 := a/4*{-1, -1, 1};
n3 := a/4*{-1, 1, -1};
n4 := a/4*{1, -1, -1};
k = {k1, k1, k1};
g0[k_] = Exp[I*(k.n1)] + Exp[I*(k.n2)] + Exp[I*(k.n3)] + Exp[I*(k.n4)];
g1[k_] = Exp[I*(k.n1)] - Exp[I*(k.n2)] - Exp[I*(k.n3)] + Exp[I*(k.n4)];
g2[k_] = Exp[I*(k.n1)] - Exp[I*(k.n2)] + Exp[I*(k.n3)] - Exp[I*(k.n4)];
g3[k_] = Exp[I*(k.n1)] + Exp[I*(k.n2)] - Exp[I*(k.n3)] - Exp[I*(k.n4)];
v0 = vss\[Sigma];
v1 = 1/Sqrt[3]*vsp\[Sigma];
v2 = 1/3*vpp\[Sigma] - 2/3*vpp\[Pi];
v3 = 1/3*vpp\[Sigma] + 1/3*vpp\[Pi];

H [k_] = {{esa, 0, 0, 0, -v0*g0[k], v1*g1[k], v1*g2[k], v1*g3[k]},
        {0, epa, 0, 0, -v1*g1[k], v2*g0[k], v3*g3[k], v3*g2[k]},
         {0, 0, epa, 0, -v1*g2[k], v3*g3[k], v2*g0[k], v3*g1[k]},
         {0, 0, 0, epa, -v1*g3[k], v3*g2[k], v3*g1[k], v2*g0[k]},
   {-v0*g0[k]\[Conjugate], -v1*g1[k]\[Conjugate], -v1*
     g2[k]\[Conjugate], -v1*g3[k]\[Conjugate], esb, 0, 0, 0},
   {v1*g1[k]\[Conjugate], v2*g0[k]\[Conjugate], v3*g3[k]\[Conjugate], 
    v3*g2[k]\[Conjugate], 0, epb, 0, 0},
   {v1*g2[k]\[Conjugate], v3*g3[k]\[Conjugate], v2*g0[k]\[Conjugate], 
    v3*g1[k]\[Conjugate], 0, 0, epb, 0},
   {v1*g3[k]\[Conjugate], v3*g2[k]\[Conjugate], v3*g1[k]\[Conjugate], 
    v2*g0[k]\[Conjugate], 0, 0, 0, epb}};

Now, I want that k1 should be real, so that the eigenvalues of H are real. How do I put the condition that k1 is real, and therefore plot the 8 eigenvalues as a function of k1.