Trouver l'équation d'une surface à partir de données

**mitch4444** · 25/03/2013, 12h40

Bonjour,

Je suis dans une impasse.
J'essaye en vain de trouver l'équation d'une surface que j'ai dessiné dans Excel dont la base de données est celle-ci :

	-90	-80	-70	-60	-50	-40	-30	-20	-10	0
0	8170	8170	8170	8170	8170	8170	8170	8170	8170	8170
10	8100	8220	8340	8440	8540	8620	8680	8730	8760	8770
20	7950	8180	8390	8590	8760	8900	9020	9110	9160	9170
30	7750	8050	8340	8600	8830	9020	9180	9300	9360	9390
40	7450	7820	8160	8470	8740	8980	9160	9300	9380	9400
50	7060	7470	7840	8190	8490	8750	8960	9110	9190	9220
60	6570	7000	7390	7760	8080	8350	8570	8730	8820	8840
70	5980	6410	6820	7180	7510	7790	8010	8160	8240	8270
80	5320	5740	6130	6490	6800	7060	7270	7400	7480	7500
90	4600	4980	5350	5670	5950	6190	6360	6470	6530	6550

A gauche (de 0 à 90) la variable x et en haut (de -90 à 0) la variable y. Pour tout ces couples de points (x,y) j'ai une valeur z.

Cela ressemble à une équation polynomiale du deuxième degré.

Avez-vous une idée à me proposer?

Merci d'avance!

**FLB** · 25/03/2013, 13h32

Salut,
Je te conseille de regarder la méthode polyfitn, disponible sur le FEX, qui te permettra de trouver les coefficients d'interpolation de ton polynôme.

**le fab** · 25/03/2013, 13h43

Salut

la fonction appropriée dans matlab est lsqcurvefit.
tu va avoir besoin également d'utiliser la fonction meshgrid pour préparer correctement tes données
et le plus important et de savoir (supposer) la forme de la fonction souhaité: quand tu parles équation polynomiale du deuxième degré, je suppose que tu fais référence à une forme quadratique.
fabien

**Jerome Briot** · 25/03/2013, 13h54

Voir aussi ici : Ajustement aux moindres carrés de courbes et de surfaces

**mitch4444** · 25/03/2013, 14h17

Envoyé par FLB

Salut,
Je te conseille de regarder la méthode polyfitn, disponible sur le FEX, qui te permettra de trouver les coefficients d'interpolation de ton polynôme.

Tout d'abord merci pour ta réponse.

J'avais déjà passé une heure à essayer de le faire fonctionner dans mon cas mais je trouve pas du tout la procédure à suivre.
Est-ce que tu aurais une piste à me donner à partir des données que j'ai cité ci-dessus?

Merci

**mitch4444** · 25/03/2013, 14h20

J'imagine que vos trois méthodes sont possibles mais j'ai vraiment du mal à les concrétiser..
Quelle ligne dois-je noter?
J'ai déjà réussi à faire apparaitre la surface à l'aide de meshgrid et surf mais après...

Merci en tout cas pour vos réponses

**le fab** · 25/03/2013, 14h47

salut

bon, déjà si tu as dessiner ta surface et utiliser meshgrid, tu as déjà fais une bonne partie du travail !
la solution proposée par FLB est une fonction faite pour ce genre de problème et normalement devrait être la plus simple à mettre en oeuvre (c'est d'ailleurs probablement une encapsulation de la méthode que te conseille Dut)
sinon pour rester dans du matlab pure, prend plutot la méthode que te conseilles Dut, elle est en effet rapide et efficace (je sais pas pourquoi, je n'y pense jamais quand j'ai un pb comme ça à résoudre, je me tourne toujours ves les fonctions d'optimisation)

lit bien le tutoriel de ol9245, tout est très bien expliqué
ça parait pas simple au début, mais ça se fait bien

après juste une remarque quand même : si tu veux mettre des contraintes min et max sur tes coefficients, les méthodes polyfitn et \ ne marcheront pas et tu devra te retourner vers des fonctions d’optimisation (lsqcurvefit voir fmincon pour intégrer plus de contraintes sur les coefficients)

**Jerome Briot** · 25/03/2013, 14h56

Envoyé par le fab

lit bien le tutoriel de ol9245, tout est très bien expliqué
ça parait pas simple au début, mais ça se fait bien

Voici un exemple :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
y = -90:10:0;
x = 0:10:90;
 
[x,y] = meshgrid(y,x);
 
z = [8170 8170 8170 8170 8170 8170 8170 8170 8170 8170
    8100 8220 8340 8440 8540 8620 8680 8730 8760 8770
    7950 8180 8390 8590 8760 8900 9020 9110 9160 9170
    7750 8050 8340 8600 8830 9020 9180 9300 9360 9390
    7450 7820 8160 8470 8740 8980 9160 9300 9380 9400
    7060 7470 7840 8190 8490 8750 8960 9110 9190 9220
    6570 7000 7390 7760 8080 8350 8570 8730 8820 8840
    5980 6410 6820 7180 7510 7790 8010 8160 8240 8270
    5320 5740 6130 6490 6800 7060 7270 7400 7480 7500
    4600 4980 5350 5670 5950 6190 6360 6470 6530 6550];
 
M = [ones(numel(z),1) x(:) y(:) x(:).^2 x(:).*y(:) y(:).^2];
 
K = M\z(:);
 
zi = zeros(size(z));
zi(:) = M*K;
 
figure
plot3(x,y,z,'r*')
hold on
surf(-90:10:0,0:10:90,zi);

**mitch4444** · 25/03/2013, 15h16

Je tente donc la méthode de l'ajustement aux moindres carrés.

Je pose X=0:10:90 et Y=-90:10:0 et j'importe mes données de Z qui donne une matrice :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
8170	8170	8170	8170	8170	8170	8170	8170	8170	8170
8100	8220	8340	8440	8540	8620	8680	8730	8760	8770
7950	8180	8390	8590	8760	8900	9020	9110	9160	9170
7750	8050	8340	8600	8830	9020	9180	9300	9360	9390
7450	7820	8160	8470	8740	8980	9160	9300	9380	9400
7060	7470	7840	8190	8490	8750	8960	9110	9190	9220
6570	7000	7390	7760	8080	8350	8570	8730	8820	8840
5980	6410	6820	7180	7510	7790	8010	8160	8240	8270
5320	5740	6130	6490	6800	7060	7270	7400	7480	7500
4600	4980	5350	5670	5950	6190	6360	6470	6530	6550

Ensuite je copie colle :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
X = X(:), Y=Y(:), Z=Z(:) ;
nx = length(X) ; 
U = ones(nx, 1) ;
M = [U X Y X.^2 X.*Y Y.^2] ;

J'obtiens un M =

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
M =
 
       1            -90              0           8100              0              0       
       1            -80             10           6400           -800            100       
       1            -70             20           4900          -1400            400       
       1            -60             30           3600          -1800            900       
       1            -50             40           2500          -2000           1600       
       1            -40             50           1600          -2000           2500       
       1            -30             60            900          -1800           3600       
       1            -20             70            400          -1400           4900       
       1            -10             80            100           -800           6400       
       1              0             90              0              0           8100

Il suffit donc plus qu'à faire le calcul :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

>> K = M\Z

Mais voici une erreur qui apparait :

??? Error using ==> mldivide
Matrix dimensions must agree.

À partir de ce stade, je ne comprends plus trop ce que je fais. Apparemment il y a une erreur car les deux matrices n'ont pas les mêmes dimensions. Comment arranger le problème?