Discussion :

[Taiby]

Bonjour,

dans un programme j'utilise la fonction Griddata pour obtenir les points interpolés entre mes points de mesure.
Le petit soucis est que ma carte de points de mesure n'est pas forcément pleine, je veux dire que sur une couronne de points il peut y a voir une zone où il n'y a pas de points de mesure. La fonction Griddata ne fait pas la différence et fait l’interpolation même aux endroits où il n'y a pas de points. Lorsque je trace la combinaison points de mesure et points interpolés sur la même figure il y a fatalement des zones où il n'y a pas de points de mesure et pourtant il y a les points interpolés (voir extrait de figure joint).



Comment puis-je "forcer" la fonction griddata à se limiter dans son interpolation ?
Dois-je passer par une autre fonction équivalent ou peut-être faire un post-traitement pour "effacer" les zones non désirées ?

Voici la partie du code où j'utilise la fonction Griddata :

Code :

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

1
2
3
4
5
6
[zzz,i]=unique(a+sqrt(-1)*b);  %suppression de tous les points (x,y) en doublon
a=a(i);                        %pour avoir une interpolation sans warning.
b=b(i);
c=c(i);
 
zi = griddata(a,b,c,xi,yi,'v4');

Merci pour votre aide.

[le fab]

salut

pour résoudre un pb équivalent, j'ai :
- délimité un contour convex de mes points réel via la fonction convhull
- puis découpé le résultat en enlevant les points en dehors de ce contour via la fonction inpolygon

c'est pas parfait (à cause du contour convexe qui n'épouse pas tes points au plus près) mais c'est déjà pas mal

Fab

[Taiby]

Merci pour ta réponse Fab,

j'ai appliqué la fonction convhull à mes données et les je les ai tracées mais l'enveloppe résultante est externe à tous mes points alors que dans mes traitements il s'agit de portions de couronnes. Il faudrait que le contour des points soit total et rentre au centre de ma couronne.

Existe t-il d'ailleurs une fonction de contour ?

[Jerome Briot]

Pourrais-tu nous fournir un jeu de données dans un fichier .mat (lui-même dans une archive zip ou rar) pour pouvoir faire quelques tests ?

[le fab]

Il faudrait que le contour des points soit total et rentre au centre de ma couronne.

donc il te faut un contour non convexe
pas de fonction toutes faites à ma connaissance ...

[Taiby]

Il me faut en effet un mixage concave/convexe... J'ai bien peur que ma requête soit biscornue


J'ai réalisé un extrait de mon code juste pour permettre de tracer un secteur représentatif de la carte finale, je l'insère dans ce message en espérant que vous pourrez m'aider.

[Jerome Briot]

Tu peux obtenir un maillage concave en utilisant l'algorithme Alpha Shapes : Alpha shapes by Jonas Lundgren

Ce qui donne concrètement avec MATLAB :

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
28
29
30
31
32
%chargement des données
a = [0;0;0;0;3.14;3.84;4.54;5.23;5.93;7.84;9.59;11.33;13.07;14.82;15.63;19.1;22.57;26.05;28.47;29.52;33.65;37.62;38.82;44;44.46;51.3;55;58.14;65;75;83.56;85;96.42;109.27;114.91;130.23];
b = [-150;-130;-110;-90;-89.95;-109.93;-129.92;-149.91;-169.9;-89.66;-109.58;-129.51;-149.43;-169.35;-88.63;-108.33;-128.03;-147.72;-106.25;-167.42;-125.57;-103.37;-144.89;-164.21;-122.16;-140.95;-95.26;-159.75;-112.58;-129.9;-99.59;-147.22;-114.91;-130.23;-96.42;-109.27];
 
%chargement des paramètres:
pGmin=0.5;
pGmax=50;
pNiso=20;
pNray=10;%17;
pNang=30;%141;
pRmin=90;
pRmax=170;
pPang=270;
pDang=320;
 
%calcul des maillages:
tta = linspace(pPang*pi/180,pDang*pi/180,pNang)'; 
r = linspace(pRmin,pRmax,pNray)';
[ri,ttai] = meshgrid(r,tta);
xi = ri.* cos (ttai);
yi = ri.* sin (ttai);
 
[V,S] = alphavol([a b],30);
 
figure
 
patch('vertices',[a b],'faces',S.tri,'facevertexcdata',S.vol, ...
    'cdatamapping','scaled','facecolor','flat');
hold on
plot(a,b,'r*',xi(:),yi(:),'bo')
 
axis equal

Voir la première image attachée.

Il "suffit" ensuite tester la position des couples (xi,yi) par rapport aux triangles générés avec alphavol pour ne conserver que ceux situés à l'intérieur du domaine.

Une autre idée qui me vient en tête pour simplifier la détection/
Pourquoi ne par "déformer" ton domaine pour le ramener sur une grille rectangulaire ?

Comme ceci :

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
28
29
30
31
32
33
34
35
36
%chargement des données
a = [0;0;0;0;3.14;3.84;4.54;5.23;5.93;7.84;9.59;11.33;13.07;14.82;15.63;19.1;22.57;26.05;28.47;29.52;33.65;37.62;38.82;44;44.46;51.3;55;58.14;65;75;83.56;85;96.42;109.27;114.91;130.23];
b = [-150;-130;-110;-90;-89.95;-109.93;-129.92;-149.91;-169.9;-89.66;-109.58;-129.51;-149.43;-169.35;-88.63;-108.33;-128.03;-147.72;-106.25;-167.42;-125.57;-103.37;-144.89;-164.21;-122.16;-140.95;-95.26;-159.75;-112.58;-129.9;-99.59;-147.22;-114.91;-130.23;-96.42;-109.27];
 
%chargement des paramètres:
pGmin=0.5;
pGmax=50;
pNiso=20;
pNray=10;%17;
pNang=30;%141;
pRmin=90;
pRmax=170;
pPang=270;
pDang=320;
 
%calcul des maillages:
tta = linspace(pPang*pi/180,pDang*pi/180,pNang)'; 
r = linspace(pRmin,pRmax,pNray)';
[ri,ttai] = meshgrid(r,tta);
xi = ri.* cos (ttai);
yi = ri.* sin (ttai);
 
[th,r] = cart2pol(a,b);
 
[thi,ri] = cart2pol(xi(:),yi(:));
 
[V,S] = alphavol([th r/100],5);
 
figure
 
patch('vertices',[th r/100],'faces',S.tri,'facevertexcdata',S.vol, ...
    'cdatamapping','scaled','facecolor','flat');
hold on
plot(th,r/100,'r*',thi(:),ri(:)/100,'bo')
 
axis equal

Voir deuxième image attachée.

Cela devrait permettre, pour chaque triangle, de ne tester que les points qui ton compris dans le rectangle l'englobant et donc d'accélérer le processus.
Il faut bien comprendre que, une fois générée, la table de connectivité du maillage est indépendante de la position des nœuds.

A tester...