Skip to content

Instantly share code, notes, and snippets.

Last active December 21, 2015 14:59
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
Star You must be signed in to star a gist
What would you like to do?
visual_dr creates and prints a web report for a dual regression run with FSL. The figures print exactly the size of the open figure 1, so you might want to open and size it to your liking before running the script. visual_dr is the main script to run, and it uses print_dr_report to create a web report, and imoverlay to overlay maps on the backgr…
function print_dr_report(imgfolder)
% This function takes a folder of images produced by visual_dr
% and prints and html report page
filey = fopen([imgfolder 'index.html' ],'w');
% Print top of html file
fprintf(filey,'%s\n','<html><body><h2>Dual Regression Report</h2>');
% Get all images in folder
images = dir([ imgfolder '*.png' ]);
for i=1:length(images)
img = images(i).name;
fprintf(filey,'%s\n',[ '<img src="' img '" \>' ]);
function visual_dr(DRdir,bgimg)
% This function will create a web report to visualize dual regression
% results. The figures print exactly the size of the open figure 1,
% so you might want to open and size it to your liking before running
% the script.
% bgimg is the background image (must be same size as dr output)
% DRdir is the dual regression directory
% Will be web report in the dual regression directory
% REQUIREMENTS - imoverlay:
% Read in the background image
bg = spm_read_vols(spm_vol(bgimg));
bgnorm = zeros(size(bg));
% Set output image directory
outdir = [DRdir '/report/' ];
for b=1:size(bg,1)
bgnorm(b,:,:) = (bg(b,:,:) - mean(bg(:))) / std(bg(:));
% Convert to between 0 and 1
bgnorm = (bgnorm - min(bgnorm(:))) ./ (max(bgnorm(:) - min(bgnorm(:))));
% Read in the number of contrasts and names
filey = fopen([ DRdir '/design.con' ],'r');
contrasts = textscan(filey,'%s','delimiter','\n');
numCons = 0;
for c=1:length(contrasts{1})
if regexp(contrasts{1}{c},'/ContrastName*') == 1
numCons = numCons + 1;
conNameIdx = regexp(contrasts{1}{c},'"*"');
cons{numCons} = contrasts{1}{c}(conNameIdx(1)+1:conNameIdx(2)-1);
% Get the number of corrected p contrast files
files = dir([ DRdir '/dr_stage3_ic00*_tfce_corrp_tstat1.nii.gz']);
% For each thresholded map, display and slice on the background image, for
% each contrast
for m=1:length(files)
for c=1:length(cons)
fig = figure(1);
% Unzip and read in the component/contrast image
mrname = regexprep(files(m).name,'tstat1',['tstat' num2str(c) ]);
corrmap = gunzip([ DRdir '/' mrname ]);
corrmap = spm_read_vols(spm_vol([ corrmap{1} ]));
% Get the cluster center of the activation
xyz = regionprops(uint8(corrmap), 'Centroid');
if ~isempty(xyz)
if length(xyz.Centroid) == 3
xyz = floor(xyz.Centroid);
% for each xyz dimension
color = rand(3,1);
subplot(1,3,1); % z
disp = imoverlay(rot90(bgnorm(:,:,xyz(3))),rot90(corrmap(:,:,xyz(3))),color);
imagesc(disp); axis off;
title([ 'Comp #' num2str(m) ' Contrast: ' cons{c}]);
subplot(1,3,2); % x
disp = imoverlay(rot90(squeeze(bgnorm(xyz(1),:,:))),rot90(squeeze(corrmap(xyz(1),:,:))),color);
imagesc(disp); axis off
subplot(1,3,3); % y
disp = imoverlay(rot90(squeeze(bgnorm(:,xyz(2),:))),rot90(squeeze(corrmap(:,xyz(2),:))),color);
imagesc(disp); axis off
set(gcf, 'PaperPositionMode','auto')
saveas(fig,[outdir 'comp' num2str(m) '_tstat' num2str(c) '.png' ])
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment