Skip to content
Snippets Groups Projects
Registration.m 14.9 KiB
Newer Older
  • Learn to ignore specific revisions
  • bepi's avatar
    bepi committed
    %% The code is written by Behnaz Pirzamanbein, bepi@dtu.dk last version 2020.08.04
    
    bepi's avatar
    bepi committed
    %% registration
    samples = {'out_out_CD','out_out_MD','in_in_CD','in_in_MD'};
    
    load_folder = 'Results\registration\';
    save_folder = 'Results\registration\tform\';
    if ~exist(save_folder,'dir')
        mkdir(save_folder)
    end
    
    % if you want to run the measure separately you can use the for loop
    %for i = 1:4
        sample = samples{i};
        % loading data
        disp(sample)
        
        dg0_inner = imread([load_folder,'registration_s1_',sample,'_0.png']);
        dg45_inner = imread([load_folder,'registration_s1_',sample,'_45.png']);
        dg90_inner = imread([load_folder,'registration_s1_',sample,'_90.png']);
        dg180_inner = imread([load_folder,'registration_s1_',sample,'_180.png']);
    
        % estimate the resgitration T (using the generated function)
        Reg_dg0_inner = registerImages(dg0_inner,dg180_inner,30);
        Reg_dg45_inner = registerImages(dg45_inner,dg180_inner,30);
        Reg_dg90_inner = registerImages(dg90_inner,dg180_inner,30);
       
        save([save_folder,'reg_tfrom_inner_',sample],'*_inner')
    
        % loading data
        dg0_outer = imread([load_folder,'registration_s2_',sample,'_0.png']);
        dg45_outer = imread([load_folder,'registration_s2_',sample,'_45.png']);
        dg90_outer = imread([load_folder,'registration_s2_',sample,'_90.png']);
        dg180_outer = imread([load_folder,'registration_s2_',sample,'_180.png']);
    
        % estimate the resgitration T (using the generated function)
        Reg_dg0_outer = registerImages(dg0_outer,dg180_outer,30);
        Reg_dg45_outer = registerImages(dg45_outer,dg180_outer,30);
        Reg_dg90_outer = registerImages(dg90_outer,dg180_outer,30);  
    
    
        if flag_save
            save([save_folder,'reg_tfrom_outer_',sample],'*_outer')
        end
    
    bepi's avatar
    bepi committed
        
    %end
    
    %% check all against 180
    
    samples = {'out_out_CD','out_out_MD','in_in_CD','in_in_MD'}; 
    degree = {'0','45','90','180'};
    
    load_folder = 'Results\registration\tform\';
    save_folder = load_folder;
    if ~exist(save_folder,'dir')
        mkdir(save_folder)
    end
    
    % if you want to run the measure separately you can use the for loop
    %for i = 1:4
        sample = samples{i};
        inner = load([load_folder,'reg_tfrom_inner_',sample]);
        outer = load([load_folder,'reg_tfrom_outer_',sample]);
    
        MOVING = inner.dg0_inner;%imrotate(layers(:,:,1),180);
        movingRefObj = imref2d(size(inner.dg0_inner));
    
        fixedRefObj = inner.Reg_dg0_inner.SpatialRefObj;
        tform = inner.Reg_dg0_inner.Transformation;
        inner_0dg = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
        MOVING = inner.dg45_inner;%imrotate(layers(:,:,1),180);
        movingRefObj = imref2d(size(inner.dg45_inner));
    
        fixedRefObj = inner.Reg_dg45_inner.SpatialRefObj;
        tform = inner.Reg_dg45_inner.Transformation;
        inner_45dg = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
        MOVING = inner.dg90_inner;%imrotate(layers(:,:,1),180);
        movingRefObj = imref2d(size(inner.dg90_inner));
    
        fixedRefObj = inner.Reg_dg90_inner.SpatialRefObj;
        tform = inner.Reg_dg90_inner.Transformation;
        inner_90dg = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
        figure,
        subplot(221)
        imshowpair(inner_0dg,inner.dg180_inner),title('0 to 180')
        subplot(222)
        imshowpair(inner_45dg,inner.dg180_inner),title('45 to 180')
        subplot(223)
        imshowpair(inner_90dg,inner.dg180_inner),title('90 to 180')
        suptitle_modified(['inner using inner T ',sample],'none')
    
        if flag_save
            saveas(gcf,[save_folder,'inner_using_inner_T_',sample,'.fig'])
            saveas(gcf,[save_folder,'inner_using_inner_T_',sample,'.png'])
        end
    
    bepi's avatar
    bepi committed
    
        % check all outer agains 180 using inner transform T
        MOVING = outer.dg0_outer;%imrotate(layers(:,:,1),180);
        movingRefObj = imref2d(size(outer.dg0_outer));
    
        fixedRefObj = inner.Reg_dg0_inner.SpatialRefObj;
        tform = inner.Reg_dg0_inner.Transformation;
        outer_0dg = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
        MOVING = outer.dg45_outer;%imrotate(layers(:,:,1),180);
        movingRefObj = imref2d(size(outer.dg45_outer));
    
        fixedRefObj = inner.Reg_dg45_inner.SpatialRefObj;
        tform = inner.Reg_dg45_inner.Transformation;
        outer_45dg = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
        MOVING = outer.dg90_outer;%imrotate(layers(:,:,1),180);
        movingRefObj = imref2d(size(outer.dg90_outer));
    
        fixedRefObj = inner.Reg_dg90_inner.SpatialRefObj;
        tform = inner.Reg_dg90_inner.Transformation;
        outer_90dg = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
        figure,
        subplot(221)
        imshowpair(outer_0dg,outer.dg180_outer),title('0 to 180')
        subplot(222)
        imshowpair(outer_45dg,outer.dg180_outer),title('45 to 180')
        subplot(223)
        imshowpair(outer_90dg,outer.dg180_outer),title('90 to 180')
        suptitle_modified(['outer using inner T ',sample],'none')
    
        if flag_save
            saveas(gcf,[save_folder,'outer_using_inner_T_',sample,'.fig'])
            saveas(gcf,[save_folder,'outer_using_inner_T_',sample,'.png'])
        end
    
    bepi's avatar
    bepi committed
    %end
    
    close all
    
    %% apply the inner T to the layers for measurement
    
    samples = {'out_out_CD','out_out_MD','in_in_CD','in_in_MD'}; 
    
    load_folder_Reg = 'Results\registration\tform\';
    load_folder_layers = 'Results\layer\';
    save_folder =  'Results\registration\layer\';
    if ~exist(save_folder,'dir')
        mkdir(save_folder)
    end
    
    % if you want to run the measure separately you can use the for loop
    %for i = 1:4
        clearvars -except samples i flag_data load_folder_Reg load_folder_layers save_folder
        close all
        
        sample = samples{i};
        fprintf('Load Data %s ...\n',sample)
        
        % for registrtaion
        inner = load([load_folder_Reg,'reg_tfrom_inner_',sample]);
        % for visualization
        outer = load([load_folder_Reg,'reg_tfrom_outer_',sample]);
    
        inner_layers_0 = load([load_folder_layers,'inner_layers_',sample,'_0']);
        outer_layers_0 = load([load_folder_layers,'outer_layers_',sample,'_0']);
    
        inner_layers_45 = load([load_folder_layers,'inner_layers_',sample,'_45']);
        outer_layers_45 = load([load_folder_layers,'outer_layers_',sample,'_45']);
    
        inner_layers_90 = load([load_folder_layers,'inner_layers_',sample,'_90']);
        outer_layers_90 = load([load_folder_layers,'outer_layers_',sample,'_90']);
    
        inner_layers_180 = load([load_folder_layers,'inner_layers_',sample,'_180']);
        outer_layers_180 = load([load_folder_layers,'outer_layers_',sample,'_180']);
    
        for s = 1:2
            fprintf('Registrtaion of Surface %d\n',s)
            
            dg0_inner = inner_layers_0.layers(:,:,s);
            dg45_inner = inner_layers_45.layers(:,:,s);
            dg90_inner = inner_layers_90.layers(:,:,s);
            dg180_inner = inner_layers_180.layers(:,:,s);
    
            dg0_outer = outer_layers_0.layers(:,:,s);
            dg45_outer = outer_layers_45.layers(:,:,s);
            dg90_outer = outer_layers_90.layers(:,:,s);
            dg180_outer = outer_layers_180.layers(:,:,s);
    
           
            MOVING = dg0_inner;
            movingRefObj = imref2d(size(dg0_inner));
    
            fixedRefObj = inner.Reg_dg0_inner.SpatialRefObj;
            tform = inner.Reg_dg0_inner.Transformation;
            [inner_0dg(:,:,s), inner_0dg_Ref] = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
            MOVING = dg45_inner;
            movingRefObj = imref2d(size(dg45_inner));
    
            fixedRefObj = inner.Reg_dg45_inner.SpatialRefObj;
            tform = inner.Reg_dg45_inner.Transformation;
            [inner_45dg(:,:,s), inner_45dg_Ref] = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
            MOVING = dg90_inner;
            movingRefObj = imref2d(size(dg90_inner));
    
            fixedRefObj = inner.Reg_dg90_inner.SpatialRefObj;
            tform = inner.Reg_dg90_inner.Transformation;
            [inner_90dg(:,:,s), inner_90dg_Ref] = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
            inner_180dg(:,:,s) = dg180_inner;
    
            MOVING = dg0_outer;
            movingRefObj = imref2d(size(dg0_outer));
    
            fixedRefObj = inner.Reg_dg0_inner.SpatialRefObj;
            tform = inner.Reg_dg0_inner.Transformation;
            [outer_0dg(:,:,s), outer_0dg_Ref] = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
            MOVING = dg45_outer;
            movingRefObj = imref2d(size(dg45_outer));
    
            fixedRefObj = inner.Reg_dg45_inner.SpatialRefObj;
            tform = inner.Reg_dg45_inner.Transformation;
            [outer_45dg(:,:,s), outer_45dg_Ref] = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
            MOVING = dg90_outer;
            movingRefObj = imref2d(size(dg90_outer));
    
            fixedRefObj = inner.Reg_dg90_inner.SpatialRefObj;
            tform = inner.Reg_dg90_inner.Transformation;
            [outer_90dg(:,:,s), outer_90dg_Ref] = imwarp(MOVING, movingRefObj, tform, 'OutputView', fixedRefObj);
    
            outer_180dg(:,:,s) = dg180_outer;
    
            figure,
            subplot(231)
            imshowpair(inner_0dg(:,:,s),inner.dg180_inner),title('0 to 180')
            ylabel('inner using inner T ','Interpreter','none')
            subplot(232)
            imshowpair(inner_45dg(:,:,s),inner.dg180_inner),title('45 to 180')
            subplot(233)
            imshowpair(inner_90dg(:,:,s),inner.dg180_inner),title('90 to 180')
    
            subplot(234)
            imshowpair(outer_0dg(:,:,s),outer.dg180_outer),title('0 to 180')
            ylabel('outer using inner T ','Interpreter','none')
            subplot(235)
            imshowpair(outer_45dg(:,:,s),outer.dg180_outer),title('45 to 180')
            subplot(236)
            imshowpair(outer_90dg(:,:,s),outer.dg180_outer),title('90 to 180')
            suptitle_modified(sample,'none')
    
            if flag_save
                saveas(gcf,sprintf('%sReg_s%d_%s.png',save_folder,s,sample))
            end
        end
        if flag_save
            disp('Save results ...')
            layers = inner_0dg;
            Ref = inner_0dg_Ref;
            save([save_folder,'Reg_inner_layers_',sample,'_0'],'layers','Ref')
            layers = inner_45dg;
            Ref = inner_45dg_Ref;
            save([save_folder,'Reg_inner_layers_',sample,'_45'],'layers','Ref')
            layers = inner_90dg;
            Ref = inner_90dg_Ref;
            save([save_folder,'Reg_inner_layers_',sample,'_90'],'layers','Ref')
            layers = inner_180dg;
            save([save_folder,'Reg_inner_layers_',sample,'_180'],'layers')
    
            layers = outer_0dg;
            Ref = outer_0dg_Ref;
            save([save_folder,'Reg_outer_layers_',sample,'_0'],'layers','Ref')
            layers = outer_45dg;
            Ref = outer_45dg_Ref;
            save([save_folder,'Reg_outer_layers_',sample,'_45'],'layers','Ref')
            layers = outer_90dg;
            Ref = outer_90dg_Ref;
            save([save_folder,'Reg_outer_layers_',sample,'_90'],'layers','Ref')
            layers = outer_180dg;
            save([save_folder,'Reg_outer_layers_',sample,'_180'],'layers')
    
    bepi's avatar
    bepi committed
        end
    %end
    
    close all 
    
    %% apply the inner T to the mesh points and return the mesh back
    
    samples = {'out_out_CD','out_out_MD','in_in_CD','in_in_MD'}; 
    
    load_folder_Reg = 'Results\registration\tform\';
    load_folder_layers = 'Results\layer\';
    save_folder = 'Results\registration\layer\';
    if ~exist(save_folder,'dir')
        mkdir(save_folder)
    end
    
    % if you want to run the measure separately you can use the for loop
    %for i = 1:4
        clearvars -except samples i flag_data load_folder_Reg load_folder_layers save_folder
        close all
        
        sample = samples{i};
        fprintf('Load Data %s ...\n',sample)
        
        % for registrtaion
        inner = load([load_folder_Reg,'reg_tfrom_inner_',sample]);
        % for visualization
        outer = load([load_folder_Reg,'reg_tfrom_outer_',sample]);
    
        inner_layers_0 = load([load_folder_layers,'inner_layers_',sample,'_0']);
        outer_layers_0 = load([load_folder_layers,'outer_layers_',sample,'_0']);
    
        inner_layers_45 = load([load_folder_layers,'inner_layers_',sample,'_45']);
        outer_layers_45 = load([load_folder_layers,'outer_layers_',sample,'_45']);
    
        inner_layers_90 = load([load_folder_layers,'inner_layers_',sample,'_90']);
        outer_layers_90 = load([load_folder_layers,'outer_layers_',sample,'_90']);
    
        inner_layers_180 = load([load_folder_layers,'inner_layers_',sample,'_180']);
        outer_layers_180 = load([load_folder_layers,'outer_layers_',sample,'_180']);
    
        for s = 1:2
            % initial transformation (rotation or flip)
            fprintf('Registrtaion of the Mesh Points Surface %d\n',s)
            
            dg0_inner = inner_layers_0.layers(:,:,s);
            dg45_inner = inner_layers_45.layers(:,:,s);
            dg90_inner = inner_layers_90.layers(:,:,s);
            dg180_inner = inner_layers_180.layers(:,:,s);
    
            dg0_outer = outer_layers_0.layers(:,:,s);
            dg45_outer = outer_layers_45.layers(:,:,s);
            dg90_outer = outer_layers_90.layers(:,:,s);
            dg180_outer = outer_layers_180.layers(:,:,s);
            
            % 0 degree
            [tmp(:,2), tmp(:,1)] = find(dg0_inner);
            tform = inner.Reg_dg0_inner.Transformation;
            [inner_0dg(:,1,s), inner_0dg(:,2,s)] = transformPointsForward(tform, tmp(:,1), tmp(:,2));
            clear tmp Reg
                    
            % 45 degree
            [tmp(:,2), tmp(:,1)] = find(dg45_inner);
            tform = inner.Reg_dg45_inner.Transformation;
            [inner_45dg(:,1,s), inner_45dg(:,2,s)] = transformPointsForward(tform, tmp(:,1), tmp(:,2));
            clear tmp Reg
            
            % 90 degree
            [tmp(:,2), tmp(:,1)] = find(dg90_inner);
            tform = inner.Reg_dg90_inner.Transformation;
            [inner_90dg(:,1,s), inner_90dg(:,2,s)] = transformPointsForward(tform, tmp(:,1), tmp(:,2));
            clear tmp Reg
            
            [inner_180dg(:,2,s), inner_180dg(:,1,s)] = find(dg180_inner);
            clear tmp Reg
    
            % 0 degree
            [tmp(:,2), tmp(:,1)] = find(dg0_outer);
            tform = inner.Reg_dg0_inner.Transformation;
            [outer_0dg(:,1,s), outer_0dg(:,2,s)] = transformPointsForward(tform, tmp(:,1), tmp(:,2));
            clear tmp Reg
            
            % 45 degree
            [tmp(:,2), tmp(:,1)] = find(dg45_outer);
            tform = inner.Reg_dg45_inner.Transformation;
            [outer_45dg(:,1,s), outer_45dg(:,2,s)] = transformPointsForward(tform, tmp(:,1), tmp(:,2));
            clear tmp Reg
    
            % 90 degree
            [tmp(:,2), tmp(:,1)] = find(dg90_outer);
            tform = inner.Reg_dg90_inner.Transformation;
            [outer_90dg(:,1,s), outer_90dg(:,2,s)] = transformPointsForward(tform, tmp(:,1), tmp(:,2));
            clear tmp Reg
            
            [outer_180dg(:,2,s), outer_180dg(:,1,s)] = find(dg180_outer);
            clear tmp Reg
    
        end
    
        if flag_save
            disp('Save results ...')
            points = inner_0dg;
            save([save_folder,'Reg_inner_points_',sample,'_0'],'points')
            points = inner_45dg;
            save([save_folder,'Reg_inner_points_',sample,'_45'],'points')
            points = inner_90dg;
            save([save_folder,'Reg_inner_points_',sample,'_90'],'points')
            points = inner_180dg;
            save([save_folder,'Reg_inner_points_',sample,'_180'],'points')
    
            points = outer_0dg;
            save([save_folder,'Reg_outer_points_',sample,'_0'],'points')
            points = outer_45dg;
            save([save_folder,'Reg_outer_points_',sample,'_45'],'points')
            points = outer_90dg;
            save([save_folder,'Reg_outer_points_',sample,'_90'],'points')
            points = outer_180dg;
            save([save_folder,'Reg_outer_points_',sample,'_180'],'points')
        end
    
    bepi's avatar
    bepi committed
    %end
    
    close all