Gs = 
        (2/((s + 4) ((s^2) + 10  s + 8) (s + 8)));

Gs = 
        (2/(256 + 416 s + 160 (s^2) + 22 (s^3) + (s^4)));
Fs = (((Subscript[K, d]) (s^2) + (Subscript[K, p]) s + (Subscript[K, i]))/s);
Gcs = 
        Together[
          ExpandAll[((Gs Fs)/(1 + Gs Fs))]];
pcs1 = Denominator[Gcs];
Print["Characteristic polynomial of the closed-loop system:", 
        pcs1];

pds = 
        ExpandAll[(s + 1/2 - 1/2 I) (s + 1/2 + 
                1/2  I)] ;
pes = (Subscript[e, 3]) (s^3) + (Subscript[e, 2 ]) (s^2) + (Subscript[e, 1]) s  +  
          (Subscript[e, 0]); 
pcs2 = ExpandAll[pds pes];
sol = 
      Solve[CoefficientList[pcs1, s] == 
          CoefficientList[pcs2, s], {(Subscript[K, p]), (Subscript[K, i]), (Subscript[e, 0]), (Subscript[e, 1]), (Subscript[e, 2]), 
          (Subscript[e, 3])}]

pefix = pes /. sol[[1]];
peshifted = ExpandAll[pefix /. s -> s - 5/2]

r1 = RouthTabulation[peshifted];
cond1 = 
        Apply[And, Map[(# > 0) &, r1]];
Reduce[cond1] // N

NFs = 
        ExpandAll[
          Numerator[Fs] /. sol[[1]]];
NFsShifted = 
        ExpandAll[NFs /. s -> s - 5/2];
r2 = RouthTabulation[NFsShifted];
cond2 = 
        Apply[And, Map[(# > 0) &, r2]];
Reduce[cond2] // N

sys = 
      DefineSystem[Gcs /. sol[[1]], LowerBounds -> {7.897}, 
        UpperBounds -> {10.0313}]