Bonjour à tous,

Je viens vers vous car travaillant actuellement sur la conception d'un imageur 3D, j'ai réalisé un algorithme permettant par le biais de la méthode itérative SART d'effectuer la reconstruction tomographique d'un objet en 3D.

Le problème est que le temps de calcul est extrêmement long (~ 5 jours pour des projection en 88*75 pixels), ceci venant du fait que la matrice de Radon (matrice permettant de passer du domaine image au domaine projection et vice versa) est de grande taille.

En utilisant le profiler, j'ai pu cibler les lignes gourmandes en temps de calcul, et c'est bien la multiplication avec la matrice de Radon qui est mise en cause.

Je me suis donc renseigné sur la parallélisation en utilisant par exemple des boucles parfor au lieu de boucles for, le problème est que j'ai vraiment du mal avec la logique de la parallélisation.

Voici la partie de mon algorithme que j'aimerai paralléliser :

De la ligne 34 (inclus) à la ligne 75 (inclus)

Code MATLAB : 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
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
%% Méthode SART
 
% Dans cette partie, nous allons utiliser une méthode de correction
% itérative de type additive.
 
 
fprintf('Méthode SART : Start\n');
 
Temporaire=zeros(L,L);
Erreur=zeros(H,ite);
pos2=1;
 
 for ligne=1:H                                                              % Balayage de l'ensemble des sinogrammes
 
    Ik1=zeros(L*L,1);                                                       % Initialisation des matrices images aux itérations k+1 et k
    Ik=zeros(L*L,1);                                                        %
    Ik2=zeros(L,L);
    Projection(:,:)=double(IMG_Sino{ligne});                                % Stockage du sinogramme = Projection_mesurée
    pos=1;
 
    for incr0=1:Nproj
        for incr=1:L
 
            Projection2(pos,1)=Projection(incr,incr0);                      %Vecteur colonne proj
 
            pos=pos+1;
        end
    end
 
 
    Projection_old=zeros(L*Nproj,1);                                        % Initialisation de la matrice projection_itération_précédente                                      
    Somme = 0;
    Err = 100;
    for  k=1:ite                                                            % Boucle définie par le nombre d'itérations choisi par l'utilisateur
 
        Soustraction = Projection2 - Projection_old;
 
        Erreur(ligne,k)=mean(Soustraction);
        Err=Erreur(ligne,k);
 
        for z=1:L*Nproj
            for a=1:L*L
 
                Somme = Somme + Matrice_Radon(z,a)*Matrice_Radon(z,a);
 
            end
 
            facteur = landa/Somme;
 
            Correction = mtimesx(facteur,Matrice_Radon_Transposee,'SPEED');                       % Formule de correction additive de l'image
 
            Correction = Correction*Soustraction;
            Ik1 = Ik + Correction;
            Ik=Ik1;                                                             % Stockage de l'image corrigée dans l'image_itération_précédente
 
             Projection_old=Matrice_Radon*Ik;                          % Calcul de la projection à partir de l'image_itération_précédente
 
        end
 
 
        pos3=1;
 
        for incr2=1:L
            for incr3=1:L
 
                Ik2(incr2,incr3)=Ik1(pos3,1);                                    % Transformation vecteur en matrice
                pos3=pos3+1;
 
            end
        end
        Coupe_SIRT(ligne,:,:)=Ik2;                                              % Stockage des coupes générées
 
    end
    fprintf('Progression : %2d/%d\n',ligne,H);                                  % Affichage de l'avancement du traitement
 end


Merci par avance, et j'espère que mon code ne vous rebutera pas trop.

Cordialement,

Alexandre.