Bonjour à tous,

Après avoir approximer le contour d'un coquillage je souhaiterais déterminer le maximum de rayon de courbure de ce dernier.
J'ai implémenté un algo mais mon souci c'est qu'il ne détecte pas le max de rayon de courbure et il me donne des valeurs aberrantes.
Par exemple pour mon graphique ci dessous man algo m'indique que le max de rayon de courbure est 1254067.537095
J'ai utilisé la formule suivante pour calculer le rayon de courbure :
RayonCourbure = ((x'*x' + y'*y')^3/2) / (x'y'' - y'x'')
Quelqu'un pourrais me dire ce qui ne vas pas dans mon algo ?

voici l'algo :
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
37
38
39
40
41
42
43
44
45
46
47
48
49
50
 
 
 
typedef struct m_Pixel
{
 
		double x;
		double y;
 
} Pixel;
 
void RayonCourbureFind( Pixel *contour, int n, double pas)
{
 
  int i,k;
  double dx,dxx,dy,dyy;
  double rayonCourbure;
  double temp;
                               /* approximation de la derivee par différence finie au point initial de la courbe */
                                dx = (contour[2].x - contour[0].x)/(pas*2.);
  				dy = (contour[2].y - contour[0].y)/(pas*2.);
 
  				dxx = (contour[2].x - 2.*contour[1].x + contour[0].x)/(pas*pas);
  				dyy = (contour[2].y - 2.*contour[1].y + contour[0].y)/(pas*pas);
 
  				/* calcul de la courbure au point initial*/
                                temp = fabs(pow(dx*dx + dy*dy, 3./2.)/(dyy*dx - dxx*dy)) ;
 
 
  for(i=1; i<n-2; ++i){
 
                                /* approximation de la derivee par différence finie  */
 
				dx = (contour[i+1].x - contour[i-1].x)/(pas*2.);
  				dy = (contour[i+1].y - contour[i-1].y)/(pas*2.);
 
  				dxx = (contour[i+1].x - 2.*contour[i].x + contour[i-1].x)/(pas*pas);
  				dyy = (contour[i+1].y - 2.*contour[i].y + contour[i-1].y)/(pas*pas);
 
				/* calcul de la courbure */
  				rayonCourbure = fabs(pow(dx*dx + dy*dy, 3./2.)/(dyy*dx - dxx*dy));                                
 
    				if(rayonCourbure > temp){
    				                     temp = rayonCourbure;
    				                     k = i;
    				                       } 	                       							    	                        
                            } 
 
   return; 
}