home *** CD-ROM | disk | FTP | other *** search
/ Celestin Apprentice 2 / Apprentice-Release2.iso / Source Code / Pascal / Libraries / GrafSys 2.0 / GrafSys 2.0 source / Transformations2.p < prev    next >
Encoding:
Text File  |  1993-06-06  |  11.7 KB  |  473 lines  |  [TEXT/PJMM]
  1. unit Transformations;
  2.  
  3. interface
  4.  
  5. {$IFC UseFixedMath = TRUE }
  6.     uses
  7.         FixedMath, Matrix;
  8.  
  9.     type
  10.  
  11.         Matrix4 = array[1..4, 1..4] of fixed;
  12.         Vector4 = array[1..4] of fixed;
  13.  
  14. {$ELSEC}
  15.     uses
  16.         Matrix;
  17.  
  18. {$ENDC}
  19.  
  20.  
  21.  
  22. (* Scale : modify transformation-Matrix T so that scaling along x,y,z is according to sx,sy,sz *)
  23.  
  24.     procedure MScale (var T: Matrix4; sx, sy, sz: real);
  25.  
  26.  
  27. (* Translate : modify tranform-Matrix T so that translation along x,y,z is according to dx,dy,dz *)
  28.  
  29.  
  30.     procedure MTranslate (var T: Matrix4; dx, dy, dz: real);
  31.  
  32. (* RotX : modify tranform-Matrix T so that  it will rotate around x-achsis according to phi *)
  33.  
  34.     procedure RotX (var T: Matrix4; phi: real);
  35.  
  36.  
  37. (* RotX : modify tranform-Matrix T so that  it will rotate around y-achsis according to phi *)
  38.  
  39.     procedure RotY (var T: Matrix4; phi: real);
  40.  
  41.  
  42. (* RotX : modify tranform-Matrix T so that  it will rotate around z-achsis according to phi *)
  43.  
  44.     procedure RotZ (var T: Matrix4; phi: real);
  45.  
  46. (* RotArbitraryAchsis : ROTATION AROUND ARBITRARY ACHSIS : *)
  47. (* p1,p2 are two Points on the Rotation-achsis (together, p1and p2 form a 3D-Line), phi is the angle *)
  48.  
  49.     procedure RotArbitraryAchsis (var T: Matrix4; p, x: Vector4; phi: real);
  50.  
  51.  
  52. (* ==================== *)
  53.  
  54. implementation
  55.  
  56. {$IFC UseFixedMath = FALSE }
  57. (* Scale : modify transformation-Matrix T so that scaling along x,y,z is according to sx,sy,sz *)
  58.  
  59.     procedure MScale (var T: Matrix4; sx, sy, sz: real);
  60.  
  61.         var
  62.             s: Matrix4;
  63.  
  64.     begin
  65.         s := Identity;
  66.         s.M[1, 1] := sx;
  67.         s.M[2, 2] := sy;
  68.         s.M[3, 3] := sz;
  69.         s.IdentityFlag := FALSE;
  70.         T := MMult(T, s);
  71.     end;
  72.  
  73.  
  74. (* Translate : modify tranform-Matrix T so that translation along x,y,z is according to dx,dy,dz *)
  75.  
  76.  
  77.     procedure MTranslate (var T: Matrix4; dx, dy, dz: real);
  78.  
  79.         var
  80.             tl: Matrix4;
  81.  
  82.     begin
  83.         tl := Identity;
  84.         tl.M[4, 1] := dx;
  85.         tl.M[4, 2] := dy;
  86.         tl.M[4, 3] := dz;
  87.         tl.IdentityFlag := FALSE;
  88.         T := MMult(T, tl);
  89.     end;
  90.  
  91.  
  92. (* RotX : modify tranform-Matrix T so that  it will rotate around x-achsis according to phi *)
  93.  
  94.     procedure RotX (var T: Matrix4; phi: real);
  95.  
  96.         var
  97.             co, si: real;
  98.             r: Matrix4;
  99.  
  100.     begin
  101.         co := cos(phi);
  102.         si := sin(phi);
  103.         r := Identity;
  104.         r.M[2, 2] := co;
  105.         r.M[2, 3] := si;
  106.         r.M[3, 2] := -si;
  107.         r.M[3, 3] := co;
  108.         r.IdentityFlag := FALSE;
  109.         T := MMult(T, r);
  110.     end;
  111.  
  112.  
  113. (* RotX : modify tranform-Matrix T so that  it will rotate around y-achsis according to phi *)
  114.  
  115.     procedure RotY (var T: Matrix4; phi: real);
  116.  
  117.         var
  118.             co, si: real;
  119.             r: Matrix4;
  120.  
  121.     begin
  122.         co := cos(phi);
  123.         si := sin(phi);
  124.         r := Identity;
  125.         r.M[1, 1] := co;
  126.         r.M[3, 1] := si;
  127.         r.M[1, 3] := -si;
  128.         r.M[3, 3] := co;
  129.         r.IdentityFlag := FALSE;
  130.         T := MMult(T, r);
  131.     end;
  132.  
  133.  
  134. (* RotX : modify tranform-Matrix T so that  it will rotate around z-achsis according to phi *)
  135.  
  136.     procedure RotZ (var T: Matrix4; phi: real);
  137.  
  138.         var
  139.             co, si: real;
  140.             r: Matrix4;
  141.  
  142.     begin
  143.         co := cos(phi);
  144.         si := sin(phi);
  145.         r := Identity;
  146.         r.M[1, 1] := co;
  147.         r.M[1, 2] := si;
  148.         r.M[2, 1] := -si;
  149.         r.M[2, 2] := co;
  150.         r.IdentityFlag := FALSE;
  151.         T := MMult(T, r);
  152.     end;
  153.  
  154. (* RotArbitraryAchsis : ROTATION AROUND ARBITRARY ACHSIS : *)
  155. (* p is a Point on the achsis, x is its Direction (together, p and x form a 3D-Line), phi is the angle *)
  156.  
  157.     procedure RotArbitraryAchsis (var T: Matrix4; p, x: Vector4; phi: real);
  158.  
  159.         var
  160.             TAA, TAAI: Matrix4;
  161.             RX, RY: Matrix4;
  162.             TRF: Matrix4;
  163.             RXT, RYT: Matrix4;
  164.             v: Vector4;
  165.             u: Vector4; (* this is the vector defined by p,x *)
  166.             alpha: real; (* rotation angle required to rotate axis to coincide with main axes*)
  167.             beta: real; (* as above *)
  168.             norm: Real; (* norm of vector x-p *)
  169.             a, b, c, d: Real; (* used to calculate sin,cos *)
  170.             cbyd, bbyd: Real; (* as above *)
  171.  
  172.     begin
  173. {step 1: transform the Achsis so that it will pass through the origin                }
  174. {            this is done by moving the first point, p, to the Origin. calculate the T-Matrix for that }
  175.  
  176.         TAA := Identity;
  177.         TAAI := Identity;
  178.         RX := Identity;
  179.         RY := Identity;
  180.  
  181.  
  182. { first, calculate the Translation Matrix to move Achsis to Rigin of coordination system }
  183.         TAA.M[4, 1] := -p[1];
  184.         TAA.M[4, 2] := -p[2];
  185.         TAA.M[4, 3] := -p[3];
  186.         TAA.IdentityFlag := FALSE;
  187.  
  188. { p[1] ist schon 1, wie TAA[4,4] }
  189. {step 2 : calculate rotation required to rotate rot-achsis to coincidence with one of the axes, here the z-achsis}
  190. {             for this we first rotate around the x-achsis to land in the z-x plane (alpa degrees) : RX }
  191. {             then, in the z-x plane, we will rotate it to coincide with z-achsis (beta degrees) : RY }
  192.  
  193.         v[1] := x[1] - p[1];
  194.         v[2] := x[2] - p[2];
  195.         v[3] := x[3] - p[3];
  196.         v[4] := 1;
  197.  
  198.         norm := sqrt(v[1] * v[1] + v[2] * v[2] + v[3] * v[3]);
  199.  
  200.         u[1] := (v[1]) / norm;
  201.         u[2] := (v[2]) / norm; (* u ist der einheitsvektor in richtung V *)
  202.         u[3] := (v[3]) / norm;
  203.         u[4] := 1;
  204.  
  205.         a := u[1];
  206.         b := u[2];
  207.         c := u[3]; (* a,b,c sind die richtungscosinuesse (?) von V  projiziert auf normal-koordinaten *)
  208.  
  209.         d := sqrt(b * b + c * c); (* sei aplpha der winkel, den V mit z-Achse bildet *)
  210.  
  211.         if d <> 0 then (* wenn d = 0, braucht nicht um x-achse gedreht zu werden *)
  212.             begin
  213.                 cbyd := c / d; (* = cos(alpha)  *)
  214.                 bbyd := b / d; (* = sin(alpha)  *)
  215.                 RX.M[2, 2] := cbyd;
  216.                 RX.M[2, 3] := bbyd;
  217.                 RX.M[3, 2] := -bbyd;
  218.                 RX.M[3, 3] := cbyd; (* durch dies wird V auf z-y-ebene gedreht *)
  219.                 RX.IdentityFlag := FALSE;
  220.             end;
  221.  
  222.         RY.M[1, 1] := d; (* d = cos(beta) *)
  223.         RY.M[1, 3] := a; (* -a = sin(beta) *)
  224.         RY.M[3, 1] := -a;
  225.         RY.M[3, 3] := d;
  226.         RY.IdentityFlag := FALSE;
  227.  
  228.         TRF := MMult(TAA, RX);
  229.         TRF := MMult(TRF, RY);
  230.  
  231. {step 3:  now we can rotate for phi degrees around z-achsis normally }
  232.  
  233.         RotZ(TRF, phi);
  234.  
  235. {step 4 : Now we have to undo all the other transformations }
  236.  
  237.         RXT := Transpose(RX); (* for rot-matrix gilt : inv(r) = transp(r) *)
  238.  
  239.         RYT := Transpose(RY);
  240.  
  241.         TAAI.M[4, 1] := p[1]; (* dies ist die inverse translate-matrix *)
  242.         TAAI.M[4, 2] := p[2];
  243.         TAAI.M[4, 3] := p[3];
  244.         TAAI.IdentityFlag := FALSE;
  245.  
  246.         TRF := MMult(TRF, RYT);
  247.         TRF := MMult(TRF, RXT);
  248.         TRF := MMult(TRF, TAAI);
  249.  
  250. {step 5 : now TRF is the matrix that contains the transformation neccessary. Update T }
  251.  
  252.         T := MMult(T, TRF);
  253.     end;
  254.  
  255.  
  256. {Now the same for Fixed-point arithmetik }
  257.  
  258. {$ELSEC}
  259. (* Scale : modify transformation-Matrix T so that scaling along x,y,z is according to sx,sy,sz *)
  260.  
  261.     procedure MScale (var T: Matrix4; sx, sy, sz: real);
  262.  
  263.         var
  264.             s: Matrix4;
  265.  
  266.     begin
  267.         s := Identity;
  268.         s.M[1, 1] := X2Fix(sx);
  269.         s.M[2, 2] := X2Fix(sy);
  270.         s.M[3, 3] := X2Fix(sz);
  271.         s.IdentityFlag := FALSE;
  272.         T := MMult(T, s);
  273.     end;
  274.  
  275.  
  276. (* Translate : modify tranform-Matrix T so that translation along x,y,z is according to dx,dy,dz *)
  277.  
  278.  
  279.     procedure MTranslate (var T: Matrix4; dx, dy, dz: real);
  280.  
  281.         var
  282.             tl: Matrix4;
  283.  
  284.     begin
  285.         tl := Identity;
  286.         tl.M[4, 1] := X2Fix(dx);
  287.         tl.M[4, 2] := X2Fix(dy);
  288.         tl.M[4, 3] := X2Fix(dz);
  289.         tl.IdentityFlag := FAlSE;
  290.         T := MMult(T, tl);
  291.     end;
  292.  
  293.  
  294. (* RotX : modify tranform-Matrix T so that  it will rotate around x-achsis according to phi *)
  295.  
  296.     procedure RotX (var T: Matrix4; phi: real);
  297.  
  298.         var
  299.             co, si: real;
  300.             r: Matrix4;
  301.  
  302.     begin
  303.         co := cos(phi);
  304.         si := sin(phi);
  305.         r := Identity;
  306.         r.M[2, 2] := X2Fix(co);
  307.         r.M[2, 3] := X2Fix(si);
  308.         r.M[3, 2] := X2Fix(-si);
  309.         r.M[3, 3] := X2Fix(co);
  310.         r.IdentityFlag := FALSE;
  311.         T := MMult(T, r);
  312.     end;
  313.  
  314.  
  315. (* RotX : modify tranform-Matrix T so that  it will rotate around y-achsis according to phi *)
  316.  
  317.     procedure RotY (var T: Matrix4; phi: real);
  318.  
  319.         var
  320.             co, si: real;
  321.             r: Matrix4;
  322.  
  323.     begin
  324.         co := cos(phi);
  325.         si := sin(phi);
  326.         r := Identity;
  327.         r.M[1, 1] := X2Fix(co);
  328.         r.M[3, 1] := X2Fix(si);
  329.         r.M[1, 3] := X2Fix(-si);
  330.         r.M[3, 3] := X2Fix(co);
  331.         r.IdentityFlag := FALSE;
  332.         T := MMult(T, r);
  333.     end;
  334.  
  335.  
  336. (* RotX : modify tranform-Matrix T so that  it will rotate around z-achsis according to phi *)
  337.  
  338.     procedure RotZ (var T: Matrix4; phi: real);
  339.  
  340.         var
  341.             co, si: real;
  342.             r: Matrix4;
  343.  
  344.     begin
  345.         co := cos(phi);
  346.         si := sin(phi);
  347.         r := Identity;
  348.         r.M[1, 1] := X2Fix(co);
  349.         r.M[1, 2] := X2Fix(si);
  350.         r.M[2, 1] := X2Fix(-si);
  351.         r.M[2, 2] := X2Fix(co);
  352.         r.IdentityFlag := FALSE;
  353.         T := MMult(T, r);
  354.     end;
  355.  
  356. (* RotArbitraryAchsis : ROTATION AROUND ARBITRARY ACHSIS : *)
  357. (* p is a Point on the achsis, x is its Direction (together, p and x form a 3D-Line), phi is the angle *)
  358.  
  359.     procedure RotArbitraryAchsis (var T: Matrix4; p, x: Vector4; phi: real);
  360.  
  361.         var
  362.             TAA, TAAI: Matrix4;
  363.             RX, RY: Matrix4;
  364.             TRF: Matrix4;
  365.             RXT, RYT: Matrix4;
  366.             v: Vector4;
  367.             u: Vector4; (* this is the vector defined by p,x *)
  368.             alpha: Fixed; (* rotation angle required to rotate axis to coincide with main axes*)
  369.             beta: Fixed; (* as above *)
  370.             norm: Fixed; (* norm of vector x-p *)
  371.             a, b, c, d: Fixed; (* used to calculate sin,cos *)
  372.             cbyd, bbyd: Fixed; (* as above *)
  373.             Rnorm: Extended;
  374.             rnx, rny, rnz: real;
  375.  
  376.     begin
  377. {step 1: transform the Achsis so that it will pass through the origin                }
  378. {            this is done by moving the first point, p, to the Origin. calculate the T-Matrix for that }
  379.  
  380.         TAA := Identity;
  381.         TAAI := Identity;
  382.         RX := Identity;
  383.         RY := Identity;
  384.  
  385.  
  386. { first, calculate the Translation Matrix to move Achsis to Rigin of coordination system }
  387.         TAA.M[4, 1] := -p[1];
  388.         TAA.M[4, 2] := -p[2];
  389.         TAA.M[4, 3] := -p[3];
  390.         TAA.IdentityFlag := FALSE;
  391.  
  392. { p[1] ist schon 1, wie TAA[4,4] }
  393. {step 2 : calculate rotation required to rotate rot-achsis to coincidence with one of the axes, here the z-achsis}
  394. {             for this we first rotate around the x-achsis to land in the z-x plane (alpa degrees) : RX }
  395. {             then, in the z-x plane, we will rotate it to coincide with z-achsis (beta degrees) : RY }
  396.  
  397.         v[1] := x[1] - p[1];
  398.         v[2] := x[2] - p[2];
  399.         v[3] := x[3] - p[3];
  400.         v[4] := FixRatio(1, 1);
  401.  
  402. (* ATTENTION: Fixed number resulted in overflow with numbers around 200. So now we *)
  403. (* use the real to calculate the norm. Note that this is not necessary in the norm calc      *)
  404. (* later on since all values have been normalized and their norm should be one                 *)
  405.  
  406.         rnx := Fix2X(v[1]);
  407.         rny := Fix2X(v[2]);
  408.         rnz := Fix2X(v[3]);
  409.         Rnorm := Sqrt(rnx * rnx + rny * rny + rnz * rnz);
  410.         norm := X2Fix(Rnorm);
  411.  
  412. (* Above was before calculated with the following statements : *)
  413. (* norm := FixMul(v[1], v[1]) + FixMul(v[2], v[2]) + FixMul(v[3], v[3]); *)
  414. (* norm := X2Fix(Sqrt(Fix2X(norm))); *)
  415. (* norm := X2Fix(sqrt(Fix2X(FixMul(v[1], v[1]) + FixMul(v[2], v[2]) + FixMul(v[3], v[3])))); *)
  416.  
  417.         u[1] := FixDiv(v[1], norm);
  418.         u[2] := FixDiv(v[2], norm); (* u ist der einheitsvektor in richtung V *)
  419.         u[3] := FixDiv(v[3], norm);
  420.         u[4] := FixRatio(1, 1);
  421.  
  422.         a := u[1];
  423.         b := u[2];
  424.         c := u[3]; (* a,b,c sind die richtungscosinuesse (?) von V  projiziert auf normal-koordinaten *)
  425.  
  426.         d := X2Fix(sqrt(Fix2X(FixMul(b, b) + FixMul(c, c)))); (* sei alpha der winkel, den V mit z-Achse bildet *)
  427.  
  428.         if d <> 0 then (* wenn d = 0, braucht nicht um x-achse gedreht zu werden *)
  429.             begin
  430.                 cbyd := FixDiv(c, d); (* = cos(alpha)  *)
  431.                 bbyd := FixDiv(b, d); (* = sin(alpha)  *)
  432.                 RX.M[2, 2] := cbyd;
  433.                 RX.M[2, 3] := bbyd;
  434.                 RX.M[3, 2] := -bbyd;
  435.                 RX.M[3, 3] := cbyd; (* durch dies wird V auf z-y-ebene gedreht *)
  436.                 RX.IdentityFlag := FALSE;
  437.             end;
  438.  
  439.         RY.M[1, 1] := d; (* d = cos(beta) *)
  440.         RY.M[1, 3] := a; (* -a = sin(beta) *)
  441.         RY.M[3, 1] := -a;
  442.         RY.M[3, 3] := d;
  443.         RY.IdentityFlag := FALSE;
  444.  
  445.         TRF := MMult(TAA, RX);
  446.         TRF := MMult(TRF, RY);
  447.  
  448. {step 3:  now we can rotate for phi degrees around z-achsis normally }
  449.  
  450.         RotZ(TRF, phi);
  451.  
  452. {step 4 : Now we have to undo all the other transformations }
  453.  
  454.         RXT := Transpose(RX); (* for rot-matrix gilt : inv(r) = transp(r) *)
  455.  
  456.         RYT := Transpose(RY);
  457.  
  458.         TAAI.M[4, 1] := p[1]; (* dies ist die inverse translate-matrix *)
  459.         TAAI.M[4, 2] := p[2];
  460.         TAAI.M[4, 3] := p[3];
  461.         TAAI.IdentityFlag := FALSE;
  462.  
  463.         TRF := MMult(TRF, RYT);
  464.         TRF := MMult(TRF, RXT);
  465.         TRF := MMult(TRF, TAAI);
  466.  
  467. {step 5 : now TRF is the matrix that contains the transformation neccessary. Update T }
  468.  
  469.         T := MMult(T, TRF);
  470.     end;
  471. {$ENDC}
  472.  
  473. end.