Discussion :

[Grasshoper]

Bonjour à toutes et à tous,

Je chercher à réaliser une fonction rapide permettant d'extraire les indices (ou la 'position') du maximum d'une matrice carrée avec pour contraintes:
1-S'il y a un seul maximum dans la matrice: indexes du maximum
2-S'il y a plusieurs maximums égaux: prendre les indexes de celui qui est le plus proche du centre de la matrice
3-S'il y a plusieurs maximums égaux à mm distance du centre de la matrice: prendre les indexes d'un de ces max au hasard

J'ai écris une partie de la solution, mais je ne suis pas satisfait de mon travail, et je suis sur qu'il existe des solutions plus rapides et moins torturées (d'autant plus que je n'ai pas encore traité le point 3 et que le code n'est pas très lisible :-/):

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
 
function ij = find_ijmax(C)
%retourne les indices du max le plus proche du centre de la matrice carrée
%C. Dans le cas de plusieurs max égaux, retourne celui le plus proche, et
%s'il en existe plusieurs au plus proche, en retourne un aléatoirement.
 
%Find indices of max(C):
[Cded,bb] = max(C);
mmC = max(max(C));
iCded = (Cded == mmC);
 
 
if(sum(iCded) > 1)
    %Several values of max(C) occured, take the nearest
    %(center of the matrix):
    m = round(size(iCded,2)/2); %take the middle position.
    %find nearest position of m.
    inds = find(iCded == 1);
    [a,b] = min(abs(m - inds)); %nearest first index from the middle position.
    ij(2) = inds(b);
 
    inds = find(C(:,ij(2)) == mmC);
    [a,b] = min(abs(m - inds));
    ij(1) = inds(b);
else
    [C,ij(2)] = max(Cded);
    ij(1) = bb(ij(2)); %position of C.
end

La matrice peut etre grande, et l'algo peut etre lancé un grand nombre de fois sur des matrices différentes, ce qui nécessite un traitement assez rapide.

Merci par avance pour vos commentaires,

Grass

[Jerome Briot]

Une idée :

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
sC = size(C);
 
[n,idx] = histc(C(:),max(C(:)));
 
[r,c] = ind2sub(sC,find(idx));
 
if n>1 % Plusieurs max
 
    d = (r-(sC(1)+1)/2).^2+(c-(sC(2)+1)/2).^2;
 
    [n,idx] = histc(d,min(d));
 
    r = r(logical(idx));
    c = c(logical(idx));
 
    if n>1 % Plusieurs min
 
        idx = randperm(numel(r));
        r = r(idx(1));
        c = c(idx(1));
 
    end
 
end

[Grasshoper]

Merci pour ce retour.
Je fais quelques tests, et la fonction me parait plus lente, est-ce du à histc ? y a t'il une alternative à cette fonction malgré le fait qu'elle soit implémentée en MEX?
merci

[Jerome Briot]

Je ne comprends pas bien ton code.

Tu considères d'abord les colonnes, puis les lignes c'est ça ?

[Grasshoper]

Oui (c'est vrai que c'est pas très clair, mais j'y travaille voir si ca peut être plus rapide (ou non?)).

[Grasshoper]

J'ai fini la fonction et je compare les 2 approches:

1ère approche avec la distance de chebyshev:

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
 
function ij = find_ijmax(C)
%retourne les indices du max le plus proche du centre de la matrice carrée
%C. Dans le cas de plusieurs max égaux, retourne celui le plus proche, et
%s'il en existe plusieurs au plus proche, en retourne un aléatoirement.
 
%Find indices of max(C):
Mat = C == max(max(C));
 
if(sum(sum(Mat)) > 1)
    m = round(size(C,2)/2); %center of the matrix.
    %recover indexes of maxs
    idx = find(Mat(:) == 1);
    %[i,j] = ind2sub(size(C,1),idx);
    i = rem(idx-1,size(C,1))+1;
    j = (idx-i)/size(C,1) + 1;
 
    %Chessboard Chebyshev distance to m:
    Cdis = max(abs(m-i),abs(m-j));
 
    %Check for min distances values:
    idx = min(Cdis) == Cdis;
    %gather indexes:
    ijs = [i(idx), j(idx)];
 
    if(size(ijs,1) > 1)
        %Several Max at equal distance:
        %randomly pick one:
        ij = ijs(ceil(rand(1)*size(ijs,1)),:);
    end
else
    %[i,j] = ind2sub(size(Mat,1),find(Mat));
    %ij = [i,j];
    %fast replacement of ind2sub:
    ij(1) = rem(find(Mat)-1,size(Mat,1))+1;
    ij(2) = (find(Mat)-ij(1))/size(Mat,1) + 1;
end

2ème approche avec histc (voir message précédent) et une petite comparaison:

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
 
clear all;
n = 100000;
tic;
for i=1:n
   C = ones(100,100);
   C(50,50) = 0;
   find_ijmax(C);
end
toc;
 
clear all;
n = 100000;
tic;
for i=1:n
   C = ones(100,100);
   C(50,50) = 0;
   find_ijmax2(C);
end
toc;

La première méthode s'avère la plus rapide (43s. vs. 69s. en moyenne sur un poste standard). Je me demande s'il est possible de gagner encore un peu de temps, bien que je ne vois pas ou creuser.

[Jerome Briot]

Quelques améliorations apportées à ta solution :

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
function ij = find_ijmax(C)
%retourne les indices du max le plus proche du centre de la matrice carrée
%C. Dans le cas de plusieurs max égaux, retourne celui le plus proche, et
%s'il en existe plusieurs au plus proche, en retourne un aléatoirement.

%Find indices of max(C):
Mat = C == max(C(:));

if(sum(Mat(:)) > 1)
    m = round(size(C,2)/2); %center of the matrix.
    %recover indexes of maxs
    idx = find(Mat(:) == 1);
    %[i,j] = ind2sub(size(C,1),idx);
    i = rem(idx-1,size(C,1))+1;
    j = (idx-i)/size(C,1) + 1;
    
    %Chessboard Chebyshev distance to m:
    Cdis = max(abs(m-i),abs(m-j));

    %Check for min distances values:
    idx = min(Cdis) == Cdis;
    %gather indexes:
    ijs = [i(idx), j(idx)];
    
    if(size(ijs,1) > 1)
        %Several Max at equal distance:
        %randomly pick one:
        ij = ijs(ceil(rand(1)*size(ijs,1)),:);
    end
else
    %[i,j] = ind2sub(size(Mat,1),find(Mat));
    %ij = [i,j];
    %fast replacement of ind2sub:
    fMat = find(Mat);
    ij(1) = rem(fMat-1,size(Mat,1))+1;
    ij(2) = (fMat-ij(1))/size(Mat,1) + 1;
end

[Grasshoper]

Merci, le gain de temps est d'en gros 50%, ce qui est bien