Identifying plane within pyramid via system of equations

I need to identify the location of the plane within a pyramid. What i have is a cross-sectional image, taken from a pyramid, which was cut horizontally. I can get the 4 edge lengths, the 4 corner angles and I also know the 4 equations of the planes, which form the pyramid. I have 12 unknowns, namely the four coordinates of my four corners (x1,y1,z1,...,x4,y4,z4,) and 12 equations (4edge lengths+ 4 corner angle + 4 plane equations).

The attached picture shows a dashed plane. This is the one we consider as an example, to verify my code. The equations of the planes are given as E1,E2,E3,E4. the lengths of the edges are all 2 and angles are 90d

egrees.

Putting all these information into the system of equation we get this:

function F = root2d(x)

F(1) = x(1)-x(3);
F(2) = -x(5)+x(6);
F(3) = x(7)+x(9)-4;
F(4) = x(11)+x(12)-4;
F(5) = sqrt((x(1)-x(4))^2+(x(2)-x(5))^2+(x(3)-x(6))^2)-2; 
F(6) = sqrt((x(4)-x(7))^2+(x(5)-x(8))^2+(x(6)-x(9))^2)-2; 
F(7) = sqrt((x(7)-x(10))^2+(x(8)-x(11))^2+(x(9)-x(12))^2)-2; 
F(8) = sqrt((x(10)-x(1))^2+(x(11)-x(2))^2+(x(12)-x(3))^2)-2; 
F(9)= (  (  (x(4)-x(1))*(x(10)-x(1))+(x(5)-x(2))*(x(11)-x(2))+(x(6)-x(3))*(x(12)-x(3))   )     /    ((sqrt((x(4))-x(1))^2+(x(5)-x(2))^2+(x(6)-x(3))^2)*sqrt((x(10)-x(1))^2+(x(11)-x(2))^2+(x(12)-x(3))^2)) )-cosd(90); 
F(10)=(  (  (x(1)-x(4))*(x(7)-x(4))+(x(2)-x(5))*(x(8)-x(5))+(x(3)-x(6))*(x(9)-x(6))   )     /    ((sqrt((x(1))-x(4))^2+(x(2)-x(5))^2+(x(3)-x(6))^2)*sqrt((x(7)-x(4))^2+(x(8)-x(5))^2+(x(9)-x(6))^2)) )-cosd(90) ;
F(11)=(  (  (x(4)-x(7))*(x(10)-x(7))+(x(5)-x(8))*(x(11)-x(8))+(x(6)-x(9))*(x(12)-x(9))   )     /    ((sqrt((x(4))-x(7))^2+(x(5)-x(8))^2+(x(6)-x(9))^2)*sqrt((x(10)-x(7))^2+(x(11)-x(8))^2+(x(12)-x(9))^2)) )-cosd(90) ;
F(12)=(  (  (x(7)-x(10))*(x(1)-x(10))+(x(8)-x(11))*(x(2)-x(11))+(x(9)-x(12))*(x(3)-x(12))   )     /    ((sqrt((x(7))-x(10))^2+(x(8)-x(11))^2+(x(9)-x(12))^2)*sqrt((x(1)-x(10))^2+(x(2)-x(11))^2+(x(3)-x(12))^2)) )-cosd(90);
end

where x1-3 are the coordinates of A, x4-6 is B, x7-9 is C, x10-12 is D. Since this is a nonlinear system of equation i apply the following to evaluate

lb=[1,1,1,1,1,1,1,1,1,1,1,1];
ub=[3,3,3,3,3,3,3,3,3,3,3,3];
rng default
x0 = [1,0.9,1.1,0.8,1.2,1.3,0.7,1.11,1.12,0.99,2.1,2.0];
[x,res] = lsqnonlin(@root2d,x0,lb,ub)

The result does not give me correct points. How can i make sure that the algorithm gives me the correct values, namely the coordinates of the points A,B,C and D?