comment extraire les coefficients d'une matrice 3x3

Question

Bonjour à tous, 

On posséde 6 séries (X,Y,Z,A,B,C) de N nombres 
(N>3 peut attiendre plusieurs dizaines) 
Exemple XYZ,ABC: array of array[0..2] of double 
On sait que ces séries sont reliées entre elles par la relation : 
(A,B,C) = M(X,Y,Z) 
où M est une matrice 3x3 

Comment extraire les 9 coefficients de cette matrice ? 

Matheux à vos claviers
jpg

cs_Kenavo · Answer

Salut,
Quand tu dis "ces séries sont reliées entre elles par la relation : (A,B,C) = M(X,Y,Z)", tu dis que chaque triplet d'indice n du tableau ABC est la solution du produit matriciel de la matrice M et du triplet d'indice n du tableau XYZ ? 
Parce que ça fait un drôle de mélange entre notation "mathématique" et pascal

Ken@vo
____________________
Code, Code, Codec !

colorid · Answer

C'est tout à fait cela.
En fait les valeurs A,B et C sont les coordonnées d'un point dans l'espace ABC correspondantes aux coordonnées X,Y,Z de ce point dans un espace XYZ sachant que ces deux espaces sont liés par une relation matricielle 3x3. C'est cette relation matricielle que je cherche.
jpg

cs_Kenavo · Answer

Bon, ben si j'ai pas tout oublié :

On va considérer que les points de chaque série sont notés :
(A1,B1,C1), (A2,B2,C2),....... et
(X1,Y1,Z1),(X2,Y2,Z2),.....
et la matice M:
((M11,M12,M13),
 (M21,M22,M23),
 (M31,M32,M33))

On a les equations:

Ai = M11.Xi + M12.Yi + M13.Zi
Bi = M21.Xi + M22.Yi + M23.Zi
Ci = M31.Xi + M32.Yi + M33.Zi

Ou si on prend 3 points

A1 = M11.X1 + M12.Y1 + M13.Z1
A2 = M11.X2 + M12.Y2 + M13.Z2
A3 = M11.X3 + M12.Y3 + M13.Z3

B1 = M21.X1 + M22.Y1 + M23.Z1
B2 = M21.X2 + M22.Y2 + M23.Z2
B3 = M21.X3 + M22.Y3 + M23.Z3

C1 = M31.X1 + M32.Y1 + M33.Z1
C2 = M31.X2 + M32.Y2 + M33.Z2
C3 = M31.X3 + M32.Y3 + M33.Z3

Ce qui nous donne trois systèmes de trois équations à trois inconnues.

une matrice Mxyz connue =

| X1 Y1 Z1 |
| X2 Y2 Z2 |
| X3 Y3 Z3 |

Il me semble me rapeller que si le déterminant de cette matrice est non nul, il y a une et une seule solution pour les valeurs Mij dans les systèmes d'équations ci-dessus.

Donc, à mon avis, pour résoudre ton problème il faut
1 - faire une boucle pour trouver dans la série XYZ trois points tels que la matrice
| Xi Yi Zi |
| Xj Yj Zj |
| Xk Yk Zk |
ait un déterminant non nul (fonction det)

2 - résoudre les 3 systèmes d'équations linéaires par la méthode de GAUSS  

                | Xi Yi Zi |   |Ai|
[M11,M12,M13] x | Xj Yj Zj | = |Aj|
                | Xk Yk Zk |   |Ak|

                | Xi Yi Zi |   |Bi|
[M21,M22,M23] x | Xj Yj Zj | = |Bj|
                | Xk Yk Zk |   |Bk|

                | Xi Yi Zi |   |Ci|
[M31,M32,M33] x | Xj Yj Zj | = |Cj|
                | Xk Yk Zk |   |Ck|

Ce qui donne un bout de code dans le genre :

const
  NbPoints = 3;
var
  ABC   : array[0..NbPoints-1] of vector;
  XYZ   : array[0..NbPoints-1] of vector;

procedure TForm1.Calcule;
var
  x,M1,M2,M3,MA,MB,MC : vector;
  mxy : Matrix;
  i,j,k,n : Integer;
  ft : TextFile;
begin
  for i:=0 to NbPoints-3 do
    for j:=i+1 to NbPoints-2 do
      for k:=j+1 to NbPoints-1 do
        begin
          for n:=1 to 3 do
            Mxy[1,n] := XYZ[i,n];
          for n:=1 to 3 do
            Mxy[2,n] := XYZ[j,n];
          for n:=1 to 3 do
            Mxy[3,n] := XYZ[k,n];
          if Det(Mxy,3) <> 0 then
            begin
              MA[1] := ABC[i,1];
              MA[2] := ABC[j,1];
              MA[3] := ABC[k,1];
              MB[1] := ABC[i,2];
              MB[2] := ABC[j,2];
              MB[3] := ABC[k,2];
              MC[1] := ABC[i,3];
              MC[2] := ABC[j,3];
              MC[3] := ABC[k,3];
              SystEq(Mxy,MA,M1,3);
              SystEq(Mxy,MB,M2,3);
              SystEq(Mxy,MC,M3,3);
              break;
            end;
        end;
end;


Et en plus ça a l'air de fonctionner, c'est rigolo !

Ken@vo
___________________
Code, Code, Codec !

cs_Kenavo · Answer

Euh ! J'ai juste oublié de te donner les références des sources

L'unité Algebra

Un cours sur le méthodes numériques matricielles

A+

Ken@vo
____________________
Code, Code, Codec !

colorid · Answer

Merci mille fois
Mais petit problème le lien sur "L'unité Algebra" :
http://www.sxlist.com/techref/language/delphi/swag/MATH0080.html 
ne fonctionne pas !
Aurais-tu une autre adresse pour télécharger cette unité ?
jpg

cs_Kenavo · Answer

Pas de problème pour moi. en voici des extraits : interface const MaxN=30;{You can increase it up to 100,not greater but each matrix variable would have size of sqr(MaxN)*sizeof(Real). It is possible to write unit for work with dinamically sized matrices, but i have no needs to do this. You can work with matrices with size less that MaxN, but while you work with this unit you must allocate memory for matrix MaxN x MaxN and leave rest of space unised} type vector=array[1..MaxN]of real; matrix=array[1..MaxN,1..MaxN]of real; sett=set of 1..MaxN; var algebrerr:boolean; procedure systeq(a:matrix;b:vector;var x:vector;n:integer); var i,j,k:integer; max:real; begin algebrerr:=false; { Conversion matrix to triangle } for i:=1 to n do begin max:=abs(a[i,i]);k:=i; for j:=succ(i) to n do if abs(a[j,i])>max then begin max:=abs(a[j,i]);k:=j end; if max<1E-10 then begin algebrerr:=true;exit end; if k<>i then begin for j:=i to n do begin max:=a[k,j]; a[k,j]:=a[i,j]; a[i,j]:=max; end; max:=b[k]; b[k]:=b[i]; b[i]:=max; end; for j:=succ(i) to n do a[i,j]:=a[i,j]/a[i,i]; b[i]:=b[i]/a[i,i]; for j:=succ(i) to n do begin for k:=succ(i) to n do a[j,k]:=a[j,k]-a[i,k]*a[j,i]; b[j]:=b[j]-b[i]*a[j,i]; end; end; { X calculation} x[n]:=b[n]; for i:=pred(n) downto 1 do begin max:=b[i]; for j:=succ(i) to n do max:=max-a[i,j]*x[j]; x[i]:=max; end; end; et function det(a:matrix;n:integer):real; var i,j,k:integer;d:real; begin for i:=1 to pred(n) do begin if abs(a[i,i])<1E-10 then begin det:=0.0;exit end; for j:=succ(i) to n do begin d:=a[j,i]/a[i,i]; for k:=i to n do a[j,k]:=a[j,k]-d*a[i,k]; end; end; d:=1.0; for i:=1 to n do d:=d*a[i,i]; det:=d; end; Ken@vo ____________________ Code, Code, Codec

colorid · Answer

ok Je viens d'essayer : Tout baigne !
Mais pour un seul triplet ABC et son correspondant XYZ !!!
Or je dispose de 36 valeurs pour chaque indice ce qui représente 7140 combinaisons possibles de triplets.
De surcroît ces valeurs sont des mesures physiques sujettes à une certaine dispersion
Il est bien sur possible de calculer les 9 coefficients sur les 7140 combinaisons et d'en faire la moyenne mais cette méthode ne détecte pas les possibles artéfacts...
Que proposes-tu pour extraire les coefficients optima ?
jpg

cs_Kenavo · Answer

Là, je crains que  cela ne dépasse mes compétences ! Il faudrait que tu t'adresses à des interlocuteurs compétents en mathématiques pour commencer, et ensuite éventuellement, chercher les algorithmes adaptés au traitement numérique approprié, et ensuite les adapter en Pascal s'il n'existent pas.

Moi, je m'arrête là ... dépassé !

bon courage 

Ken@vo
____________________
Code, Code, Codec !

colorid · Answer

De toutes façons tes algos me sont utiles et seront utilisés, car
il existe des méthodes dites des moindres carrés qui utilisées avec des matrices deviennent vite ingérables et ne permettent pas d'isoler les artéfacts. En fait il s'agit pour moi de trouver une bonne méthode pour trouver les triplets les plus représentatifs pour ne pas me taper 1740 calculs...
Encore merci et bon vent...
jpg

cs_jmic · Answer

Bonjour,



J'ai achet&#233; r&#233;cemment un livre sur la m&#233;thode des moindres carr&#233;s aux
&#233;ditions Herm&#233;s fait par un monsieur qui travaille &#224; l'Institut
G&#233;ographique National. (Estimation par moindres carr&#233;s, P. Sillard,
Herm&#232;s / Lavoisier)

C'st pas hyper compliqu&#233; (mais quand m&#234;me) et devrait te permettre de r&#233;soudre ton probl&#232;me au mieux

Ceci dit il doit bien exister quelque part un programme tout fait qui prend les vecteurs en entr&#233;e et sort le r&#233;sultat.

jmic

Comment extraire les coefficients d'une matrice 3x3

10 réponses

Votre réponse

Discussions similaires