Publicité
+ Répondre à la discussion
Affichage des résultats 1 à 14 sur 14
  1. #1
    Candidat au titre de Membre du Club
    Inscrit en
    août 2010
    Messages
    48
    Détails du profil
    Informations forums :
    Inscription : août 2010
    Messages : 48
    Points : 13
    Points
    13

    Par défaut Equation de plan à partir de points

    Bonjour,

    j'ai suivi ce tutoriel pour déterminer un plan a partir de plusieurs points.

    J'ai essayé de suivre au mieux les indications, mais mon resultat est peu probant :
    Code :
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    14
        X = Pts(:,1)
        Y = Pts(:,2)
        Z = Pts(:,3)
        U = ones(length(X), 1);
        M = [U X Y X.^2 X.*Y Y.^2]
        K = M\Z;
        Z;
        Z_modele = M*K;
        residus = Z-Z_modele; 
        variance_d_origine = var(Z);
        variance_expliquee = var(Z_modele);
        variance_residuelle = var(residus);
        coefficient_correlation_r2 = variance_expliquee/variance_d_origine
    Ou X contient tous les X de mes points (x1,x2,...xn), Y et Z egalement sur le meme principe.

    Je souhaiterais dessiner mon plan calculé mais je ne sais pas comment m'y prendre.

    Merci

  2. #2
    Membre Expert Avatar de davcha
    Profil pro
    Inscrit en
    avril 2004
    Messages
    1 257
    Détails du profil
    Informations personnelles :
    Âge : 33
    Localisation : France

    Informations forums :
    Inscription : avril 2004
    Messages : 1 257
    Points : 1 399
    Points
    1 399

    Par défaut

    Une autre façon de faire consiste à calculer les valeurs singulières de la matrice formée par tes points, centrés en zéro.

    En gros :

    Code :
    1
    2
    3
    4
    5
    6
    7
    M=Pts - mean(Pts)
    [u s v]=svd(M)
    n=v(:,3) % normale à l'hyperplan.
    n/=n(3)
    [X Y]=meshgrid(-10:10)
    Z=-n(1)*X-n(2)*Y
    mesh(X,Y,Z)

  3. #3
    Candidat au titre de Membre du Club
    Inscrit en
    août 2010
    Messages
    48
    Détails du profil
    Informations forums :
    Inscription : août 2010
    Messages : 48
    Points : 13
    Points
    13

    Par défaut

    Merci pour ton aide, mais il y a un probleme a la ligne : n/=n(3) et je n'ai pas suffisamment compris pour corriger moi meme

  4. #4
    Membre Expert Avatar de davcha
    Profil pro
    Inscrit en
    avril 2004
    Messages
    1 257
    Détails du profil
    Informations personnelles :
    Âge : 33
    Localisation : France

    Informations forums :
    Inscription : avril 2004
    Messages : 1 257
    Points : 1 399
    Points
    1 399

    Par défaut

    n/=n(3) c'est équivalent à n=n/n(3).

  5. #5
    Candidat au titre de Membre du Club
    Inscrit en
    août 2010
    Messages
    48
    Détails du profil
    Informations forums :
    Inscription : août 2010
    Messages : 48
    Points : 13
    Points
    13

    Par défaut

    je vois, mais j'ai une erreur de dimension en utilisant mean(). Je dois soustraire la valeur moyenne de X, Y et Z pour chaque point ?

    Par exemple : X(1) - moy(X) ?

  6. #6
    Membre Expert Avatar de davcha
    Profil pro
    Inscrit en
    avril 2004
    Messages
    1 257
    Détails du profil
    Informations personnelles :
    Âge : 33
    Localisation : France

    Informations forums :
    Inscription : avril 2004
    Messages : 1 257
    Points : 1 399
    Points
    1 399

    Par défaut

    oui.

  7. #7
    Candidat au titre de Membre du Club
    Inscrit en
    août 2010
    Messages
    48
    Détails du profil
    Informations forums :
    Inscription : août 2010
    Messages : 48
    Points : 13
    Points
    13

    Par défaut

    Merci beaucoup, cela fonctionne "presque". J'obtiens en effet un plan avec la bonne inclinaison, mais j'aurais aimé qu'ils soit proche des points, il se positionne pas toujours correctement, j'ajoute + mean(Z)

    Code :
    1
    2
     
        Z=-n(1)*X-n(2)*Y + mean(Z)

  8. #8
    Candidat au titre de Membre du Club
    Inscrit en
    août 2010
    Messages
    48
    Détails du profil
    Informations forums :
    Inscription : août 2010
    Messages : 48
    Points : 13
    Points
    13

    Par défaut

    Cette méthode donne la bonne inclinaison du plan, mais la position en Z est erronée.

    Mes modifications precedentes n'y ont rien changé

  9. #9
    Membre Expert Avatar de davcha
    Profil pro
    Inscrit en
    avril 2004
    Messages
    1 257
    Détails du profil
    Informations personnelles :
    Âge : 33
    Localisation : France

    Informations forums :
    Inscription : avril 2004
    Messages : 1 257
    Points : 1 399
    Points
    1 399

    Par défaut

    ok, une autre façon de faire...

    Code :
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    M = [Pts ones(size(Pts,1),1)];
    [u s v] = svd(M);
    v=v(:,4);
    # v(1)*x + v(2)*y + v(3)*z + v(4) = 0
    v/=v(3);
    v=[v(1:2) v(4)];
    
    # nouveaux points...
    x = [randn(1000,2) ones(1000,1)];
    x(:,3) = x*v;
    
    # tous les x sont sur ton plan.

  10. #10
    Candidat au titre de Membre du Club
    Inscrit en
    août 2010
    Messages
    48
    Détails du profil
    Informations forums :
    Inscription : août 2010
    Messages : 48
    Points : 13
    Points
    13

    Par défaut

    Au risque de jouer les rabajoie, j'obtiens une erreur de dimension a

    et la

    Merci quand meme pour ton temps

  11. #11
    Membre Expert Avatar de davcha
    Profil pro
    Inscrit en
    avril 2004
    Messages
    1 257
    Détails du profil
    Informations personnelles :
    Âge : 33
    Localisation : France

    Informations forums :
    Inscription : avril 2004
    Messages : 1 257
    Points : 1 399
    Points
    1 399

    Par défaut

    v=[v(1:2);v(4)];

    pardon.

  12. #12
    Candidat au titre de Membre du Club
    Inscrit en
    août 2010
    Messages
    48
    Détails du profil
    Informations forums :
    Inscription : août 2010
    Messages : 48
    Points : 13
    Points
    13

    Par défaut

    J'ai essayé de modifier mon code de cette maniere, mais je crois faire fausse route. Le probleme est que je n'ai qu une vague idee de ce que je poste la...

    Code :
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    M = [Pts ones(size(Pts,1),1)];
        [u s v] = svd(M);
        v=v(:,4);
        % # v(1)*x + v(2)*y + v(3)*z + v(4) = 0
        v=v/v(3);
        v=[v(1:2);v(4)];
        [X Y]=meshgrid(0:5:300);
        Z=-v(1)*X-v(2)*Y + mean(Z);
        hold on
        o = mesh(X,Y,Z,'FaceLighting','gouraud','LineWidth',0.3);
        colormap hot;

  13. #13
    Membre Expert Avatar de davcha
    Profil pro
    Inscrit en
    avril 2004
    Messages
    1 257
    Détails du profil
    Informations personnelles :
    Âge : 33
    Localisation : France

    Informations forums :
    Inscription : avril 2004
    Messages : 1 257
    Points : 1 399
    Points
    1 399

    Par défaut

    La ligne 8 est pas bonne.

    Quand tu fais la SVD de la matrice M, la 4e colonne des vecteurs propres à droite correspond à l'équation du plan, comme indiqué dans la ligne 4.

    Donc, pour obtenir Z, la 3e coordonnée, tu fais : -(v(1)*x+v(2)*y+v(4)) / v(3)
    Donc, pour simplifier le calcul, tu pars de ton vecteur v original (celui de la ligne 3), tu divises chaque composante par v(3). Ensuite tu crées un nouveau vecteur V=-[v(1:2);v(4)]. A partir de là, tout point [x y 1] multiplié par le vecteur V te donnera la 3e coordonnée, Z.

    Dans le cas des meshgrid, tu dois donc faire ceci :

    Z=-(v(1)*x+v(2)*y+v(4))/v(3)

  14. #14
    Candidat au titre de Membre du Club
    Inscrit en
    août 2010
    Messages
    48
    Détails du profil
    Informations forums :
    Inscription : août 2010
    Messages : 48
    Points : 13
    Points
    13

    Par défaut

    Merci beaucoup, j'ai enfin trouvé !

    j'avoue m'etre aidé de ce poste, mais c'est juste parce que j'ai du mal.

    Je donne quand meme la solution car il y avait un pti souci avec le vecteur v(4):

    Code :
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    M = [Pts ones(size(Pts,1),1)];
        [u s v] = svd(M);
        v=v(:,4);
        % # v(1)*x + v(2)*y + v(3)*z + v(4) = 0
        v=v/v(3);
        v=-[v(1:2);v(4)];
        [X Y]=meshgrid(0:5:300);
        Z=v(1)*X+v(2)*Y+v(3)
        hold on
        o = mesh(X,Y,Z);
    Une autre méthode (celle du lien) plus simple pour moi sans svd() :
    Code :
    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
        x = Pts(:,1)';
        y = Pts(:,2)';
        z = Pts(:,3)';
        M = [x; y; ones(1,size(Pts,1))];
        c = z/M;
        [X Y]=meshgrid(0:5:300);
        Z=c(1)*X+c(2)*Y+c(3);
        hold on
        plot3(x,y,z,'r*')
        hold on
        o = mesh(X,Y,Z)
    Conclusion : je suis vraiment bidon en maths.

    Merci beaucoup pour ton aide !

+ Répondre à la discussion
Cette discussion est résolue.

Liens sociaux

Règles de messages

  • Vous ne pouvez pas créer de nouvelles discussions
  • Vous ne pouvez pas envoyer des réponses
  • Vous ne pouvez pas envoyer des pièces jointes
  • Vous ne pouvez pas modifier vos messages
  •