Implement algorithm to convert tiff stacks to cube format

[tmbo]

Log Time

[tmbo]

Code from Florian as an orientation:

% Writes Tiff Stacks into Knossos/Oxalis Hierarchy
% Aug 2012 florian@mhlab.net
% Input: 1x1xZ tiff files
% For cubing of higher mags input images have to be already downsampled in ALL dimensions (also z!!) 
% Caches arbitrary number of images
% Deletes empty cubes after each z iteration

%% User Parameters
dirSource='C:\...';
dirTarget='C:\...';

% Infos written into Knossos.conf file
expName='TK_mag1';
mag=1;
xyzScale=[200 200 200];

% ContainerSize: Specify how many images you want to be cached at one time (maximum depends on
% your computers ram). Must be a divisor of cubeSize, e.g. 128, 64, 32, ...
containerSize=128;

% Knossos Format = 128
cubeSize=128;

%% Code
if mod(cubeSize, containerSize)~=0
    error('ERROR: ContainerSize must be a divisor of cubeSize, e.g. 128, 64, 32, 16, ...')
end
warning('off','all');

%List files of specified type in given folder
filelist=dir(dirSource);
j=0;
for i=1:size(filelist,1)
    if isempty(regexpi(filelist(i).name, '.tif'))==0
        j=j+1;
        imagelist(j,1)=filelist(i);
    end
end

% Compute Image Parameters
imgNumber=size(imagelist,1);
img=imread(fullfile(dirSource, imagelist(1).name));
imgHeight=size(img,1);
imgWidth=size(img,2);

% Compute Number of Cubes needed in every direction
cubesX=floor((imgWidth-1)/cubeSize)+1;
cubesY=floor((imgHeight-1)/cubeSize)+1;
cubesZ=floor((imgNumber-1)/cubeSize)+1;
zssRatio=cubeSize/containerSize; % Number of z substacks per cube
zssTotal=floor(imgNumber/containerSize)+1;
disp('Initializing Image Container');
imgContainer=zeros((cubesX+1)*128, (cubesY+1)*128, containerSize, 'uint8');

% Load images & write cubes
% xyz are the anchors of the cube hierarchy

for z=0:cubesZ-1
    c=0;
    zstring=strcat('z',sprintf('%04.0f',z));
    imgContainer(:,:,:)=0;

    % localSubZ iterates through the subcubes of each z-cube (number of
    % subcubes given by cubesize/containersize
    for zssL=1:zssRatio % zssL = position within local substack {1:z_ssdiv}
        zssG=((z)*zssRatio)+zssL; % zssG = global position within all substacks z {1:z_div*cubesZ)
        disp('Loading Image Container...');

        % Container is filled with images
        if zssG<zssTotal 
            for i=1:containerSize % z position within cube
                j=i+(zssL-1)*containerSize+(z*cubeSize); % absolute z position within stack
                imgPlane=imread(fullfile(dirSource, imagelist(j).name));
                imgPlane=imrotate(flipdim(imgPlane,2),90);
                imgContainer(1:imgWidth,1:imgHeight,i)=imgPlane(1:imgWidth,1:imgHeight); 
            end            
        elseif zssG==zssTotal % for the last subcube where corresponding images are available load as much as we have
            for i=1:mod(imgNumber,containerSize) % z position within cube
                j=i+(zssL-1)*containerSize+(z*cubeSize); % absolute z position within stack
                imgPlane=imread(fullfile(dirSource, imagelist(j).name)); 
                imgPlane=imrotate(flipdim(imgPlane,2),90);
                imgContainer(1:imgWidth,1:imgHeight,i)=imgPlane(1:imgWidth,1:imgHeight); 
            end            
        elseif zssG>zssTotal
            imgContainer(:,:,:)=0;            
        end              
        cuboid=zeros(cubeSize,cubeSize,containerSize,'uint8');        
        disp('Writing Cubes...');
        id=0;
        for y=0:(cubesY-1)
            for x=0:(cubesX-1)
                id=id+1;
                tic
                anchor=[(y*cubeSize)+1 (x*cubeSize)+1];
                cuboid(1:cubeSize,1:cubeSize,1:containerSize)=imgContainer(anchor(2):anchor(2)+cubeSize-1,anchor(1):anchor(1)+cubeSize-1,1:containerSize);
                xstring=strcat('x',sprintf('%04.0f',x));
                ystring=strcat('y',sprintf('%04.0f',y));                         
                cubepath=strcat(dirTarget,xstring,'\',ystring,'\',zstring,'\'); 
                cubename=strcat(expName,'_',xstring,'_',ystring,'_',zstring,'.raw');
                mkdir(cubepath);
                fid=fopen(strcat(cubepath,cubename),'a+');
                fwrite(fid, cuboid);
                fclose(fid);               
                if cuboid(:,:,:)==0
                    c=c+1;
                    emptycubes{c,1}=strcat(cubepath,cubename);
                end                
                disp(strcat('part',{' '},num2str(zssL),{' '},'of',{' '},num2str(zssRatio),{' '},'of',{' cube '}, xstring,'_',ystring,'_',zstring,' written'));
                toc
            end
        end
    end
    disp('Deleting Empty Cubes');
    [b,sx,sy]=unique(emptycubes);
    counts=accumarray(sy,1);
    cuboidcounts=[b,num2cell(counts)];
    for cb=1:size(cuboidcounts,1)
        if cuboidcounts{cb,2}==zssRatio
            delete(cell2mat(cuboidcounts(cb,1)));
            foldername=regexpi(cuboidcounts(cb,1), '(.*\\x.*\\y.*\\z.*)\\.*', 'tokens');
            rmdir(foldername{1});
        end
    end  
    if zssL==zssRatio
        progress=(z+1)/cubesZ*100;
        disp(strcat('Overall Progress ~',sprintf('% 3.0f',progress), ' %'));
    end       
end
fclose all;

% kconfig
kl_fname = strcat(dirTarget,'knossos.conf');
fid=fopen(kl_fname,'w');
fprintf(fid,'experiment name \"%s\";\n',expName);
fprintf(fid,'boundary x %d;\n',cubesX*cubeSize);
fprintf(fid,'boundary y %d;\n',cubesY*cubeSize);
fprintf(fid,'boundary z %d;\n',cubesZ*cubeSize);
fprintf(fid,'scale x %.2f;\n',xyzScale(1));
fprintf(fid,'scale y %.2f;\n',xyzScale(2));
fprintf(fid,'scale z %.2f;\n',xyzScale(3));
fprintf(fid,'magnification %d;\n',mag);
fclose(fid);
disp('Knossos.conf written');