Calcul précis de la dérivée d'une fonction pour trouver des pics

Version imprimable

Bonjour,

J'ouvre ce sujet qui est en fait un complément de mon autre sujet : lien où j'essaye de détecter les extremums d'une fonction entre deux bornes. Le problème qui se pose est le calcul de la dérivée de la fonction ou tout au moins d'une différence. Voici d'abord un code à copier/coller pour obtenir le graph :
Code:

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 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 function test A1=[23.588 29.936 27.605 26.785 24.608 27.746 25.733 26.833 28.643 28.529 26.804 26.709 24.703 24.939 32.869 28.473 28.091 27.063 26.669 32.049 28.024 28.912 30.956 24.751 26.907 29.645 27.981 32.656 30.61 25.18 29.272 25.97 28.047 27.803 24.802 25.902 29.432 26.64 28.933 20.475 -44.146 -112.87 -183.31 -253.61 -319.41 -394.98 -460.09 -528.45 -601.48 -505.61 -605.73 -683.95 -463.52 -388.24 -325.31 -250.72 -181.37 -111.72 -45.026 20.974 30.05 28.257 27.35 25.883 26.945 23.549 26.536 28.269 25.945 23.717 32.259 30.745 26.194 26.054 27.38 31.165 29.668 26.674 27.556 24.05 25.536 28.851 26.579 30.121 28.352 30.414 29.367 25.492 25.602 26.171 26.842 25.852 27.751 29.553 31.293 25.307 30.274 27.687 24.156 24.295]; figure(1) plot(A1) title('signal A50 : Bruit+Reflux+Sinus') xlabel('indices') ylabel('valeurs indices')
Comme vous allez le voir ce qui m'intéresse c'est uniquement les pics entre les indices 40 et 60 et plus précisément le maximum (qui est atteint en l'indice 50) et le minimum (en l'indice 52). Enfin si j'arrive déjà à détecter les pics il me suffira de conserver le plus grand maximum et le plus petit minimum.

En parcourant le forum j'ai essayé la solution de Dut (ici), récapitulée dans cette seconde fonction :
Code:

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 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 function test2 y=[23.588 29.936 27.605 26.785 24.608 27.746 25.733 26.833 28.643 28.529 26.804 26.709 24.703 24.939 32.869 28.473 28.091 27.063 26.669 32.049 28.024 28.912 30.956 24.751 26.907 29.645 27.981 32.656 30.61 25.18 29.272 25.97 28.047 27.803 24.802 25.902 29.432 26.64 28.933 20.475 -44.146 -112.87 -183.31 -253.61 -319.41 -394.98 -460.09 -528.45 -601.48 -505.61 -605.73 -683.95 -463.52 -388.24 -325.31 -250.72 -181.37 -111.72 -45.026 20.974 30.05 28.257 27.35 25.883 26.945 23.549 26.536 28.269 25.945 23.717 32.259 30.745 26.194 26.054 27.38 31.165 29.668 26.674 27.556 24.05 25.536 28.851 26.579 30.121 28.352 30.414 29.367 25.492 25.602 26.171 26.842 25.852 27.751 29.553 31.293 25.307 30.274 27.687 24.156 24.295]; x=1:size(y,1); sa=sign(diff([-inf y'])); sb=sign(diff([-inf y(end:-1:1)'])); sb=sb(end:-1:1); idx=(sa==1 & sb==1); figure plot(x,y,'b-',x(idx),y(idx),'r*')
Vous pouvez tester mais sinon voici le résultat :

http://i1240.photobucket.com/albums/...ps63353d42.png

Malheureusement entre 40 et 60 les deux minimum ne sont pas détectés, seulement le maximum ! Comment faire pour avoir tous les extremas entre ces deux bornes seulement ?

EDIT

Cette autre variante, toujours de Dut semble en fait fonctionner miraculeusement :
Code:

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 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 function test2 y=[23.588 29.936 27.605 26.785 24.608 27.746 25.733 26.833 28.643 28.529 26.804 26.709 24.703 24.939 32.869 28.473 28.091 27.063 26.669 32.049 28.024 28.912 30.956 24.751 26.907 29.645 27.981 32.656 30.61 25.18 29.272 25.97 28.047 27.803 24.802 25.902 29.432 26.64 28.933 20.475 -44.146 -112.87 -183.31 -253.61 -319.41 -394.98 -460.09 -528.45 -601.48 -505.61 -605.73 -683.95 -463.52 -388.24 -325.31 -250.72 -181.37 -111.72 -45.026 20.974 30.05 28.257 27.35 25.883 26.945 23.549 26.536 28.269 25.945 23.717 32.259 30.745 26.194 26.054 27.38 31.165 29.668 26.674 27.556 24.05 25.536 28.851 26.579 30.121 28.352 30.414 29.367 25.492 25.602 26.171 26.842 25.852 27.751 29.553 31.293 25.307 30.274 27.687 24.156 24.295]; figure plot(y) s=sign(diff([+inf y' +inf])); idx=(s(2:end).*s(1:end-1))==-1; hold on plot(find(idx),y(idx),'ro')
http://i1240.photobucket.com/albums/...ps4be2f501.jpg

Cependant je voudrais faire fonctionner ceci seulement entre les bornes 40 et 60.

25/07/2013, 09h34
abel413
Par exemple ceci ne marche pas :

Code:

s=sign(diff([+inf y(40:60)' +inf]));

Je n'obtiens pas les trois extrema entre 40 et 60. Vous voyez le problème ? Si je fais bien le travail sur toute la courbe j'obtiens certes comme suit les extrema :
Code:

1 2 3 indices_pics=find(idx)'; valeurs_indices_pics=y(idx); extrema=[indices_pics valeurs_indices_pics];
Mais les pics qui m'intéressent aux rang 49,50 et 52 sur la première colonne de extrema correspondent aux indices 26,27 et 28. Bref mon problème est un problème d'indicage tout bête ! J'ai besoin d'obtenir les pics mais avec les bons indices pour pouvoir les récupérer ou, alors, comme je voulais, n'obtenir que les pics entre deux bornes dès le départ.
25/07/2013, 10h07
Jerome Briot

Ajoute simplement la ligne suivante juste avant le hold on.

Code:

idx([1:39 61:end]) = false;

Merci Dut, effectivement je n'avais pas assez réfléchi. Malgré tout j'ai procédé d'une autre manière entre temps. Tout marche bien mais je me heurte à un souci (probablement hors-sujet) dont j'avais discuté avec Caro-Line dans un autre sujet mais je n'arrive toujours pas à faire comme il faut. Voici déjà le code utilisé :
Code:

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 37 38 39 40 41 42 43 44 45 %Visualisation de la courbe à étudier. figure(3) plot(y) xlabel('indices') ylabel('valeurs indices') s=sign(diff([+inf y' +inf])); idx=(s(2:end).*s(1:end-1))==-1; %Indices et valeurs associées des différents extrema trouvés. A ce stade on %a que les extremas sans savoir si ce sont des maxima ou des minima. On va %se restreintre ici aux extremas d'indices 45 à 55 environ. indices_pics=find(idx)'; valeurs_indices_pics=y(idx); extrema=[indices_pics valeurs_indices_pics]; extrema_retenu=find(extrema(:,1)>45 & extrema(:,1)<55 ); extrema=[indices_pics(extrema_retenu) valeurs_indices_pics(extrema_retenu)]; %Détermination des deux pics à conserver : on doit conserver le plus grand %maximum et le plus petit minimum. [plus_grand_max i_plus_grand_max]=max(extrema(:,2)); [plus_petit_min i_plus_petit_min]=min(extrema(:,2)); %Supression des lignes d'indices différents de i_plus_grand_max et %i_plus_petit_min. for i=1:size(extrema,1) if (i~=i_plus_grand_max & i~=i_plus_petit_min) extrema(i,:)=0; else continue; end end Z1=find(extrema(:,1)==0); extrema(Z1,:)=[]; %On affiche les deux extremas retenus. hold on plot(extrema(:,1),extrema(:,2),'ro')
J'ai décidé de stocker dans la matrice 2 colonnes extrema l'ensemble des extrema trouvés puis dans extrema_retenu je ne conserve que les extremas entre les bornes 45 à 55 et je remets à jour extrema. Avec l'exemple proposé plus haut, il me reste donc bien trois pics, voire 4 selon la situation. Il me reste alors à trouver le plus grand maximum et le plus petit minimum ce qui est fait avec min et max. Et voilà, ensuite je supprime à l'aide d'une boucle for les lignes de extrema dont les indices sont différents des indices min et max trouvés. Je fais quelque chose d'archaïque avec une boucle for et un if dedans : je mets à 0 les lignes qui ne m'intéressent pas puis en dehors de la boucle je recherche les 0 sur la première colonne de extrema et j'enlève les lignes correspondante. Bref cela marche mais je suis sûr qu'il y a plus simple avec setdiff comme m'en avait parlé Caro-Line mais je n'y arrive pas. En résumé je veux donc dans extrema supprimer toutes lignes d'indices autres que celles d'indices i_plus_grand_max et i_plus_grand_min.