function [ output_args ] = frrcf(parent_directory,alg,scale,threshold,iterations,steps,reduc,frequencyX,frequencyY,frequencyZ  )

%=============================Get User Input===============================

T_N_Images = length(dir([parent_directory,'\*Unmixed_Turquoise*.tif']));
Number_zStacks = length(dir([parent_directory,'\*01_Unmixed_Turquoise*.tif']));
if Number_zStacks == 0
    Number_zStacks = length(dir([parent_directory,'\*1_Unmixed_Turquoise*.tif']));
end
NumberImages = T_N_Images./Number_zStacks;
% we used Guassian smoothing
if  strcmp(alg,'Gaussian')
    smoothing = 1;
end

if strcmp(alg,'AnisoDiff_2D')
    smoothing = 2;
end

if strcmp(alg,'AnisoDiff_3D')
    smoothing = 3;
end

if strcmp(alg,'fast_AnisoDiff_3D')
    smoothing = 4;
end

%Create a directory to save resliced data
resliced_directory=strcat(parent_directory,'\Resliced RAWandLUTFRET-',num2str(threshold),'threshold');
mkdir(resliced_directory);

fid=fopen(strcat(resliced_directory,'\smoothing_settings.txt'),'w');
fprintf(fid,['scale = ',num2str(scale)]);
fprintf(fid,'\r\n');
fprintf(fid,['frequencyX = ',num2str(frequencyX)]);
fprintf(fid,'\r\n');
fprintf(fid,['frequencyY = ',num2str(frequencyY)]);
fprintf(fid,'\r\n');
fprintf(fid,['frequencyZ = ',num2str(frequencyZ)]);
fprintf(fid,'\r\n');
fprintf(fid,['Smoothing Type: ',alg]);

if strcmp(alg,'Gaussian') == 0
fprintf(fid,'\r\n');
fprintf(fid,['iterations = ',num2str(iterations)]);
fprintf(fid,'\r\n');
fprintf(fid,['steps = ',num2str(steps)]);
fprintf(fid,'\r\n');
fprintf(fid,['reduction strength = ',num2str(reduc)]);
end

%==========================Load Original Images============================

 %search and load function   function[ image3DOut,mImageOut,nImageOut ]=searchAndLoad(parent_directory,NumberImages,Number_zStacks,keyword)
 function[ image3DOut,mImageOut,nImageOut ]=searchandload(parent_directory,NumberImages,timepoints,keyword)      
 fullSearchName = strcat(parent_directory,'\*',keyword,'.tif');
        listing = dir(fullSearchName); %retrieves information on all specified files

        InfoImage=imfinfo(strcat(parent_directory,'\',listing(1).name));
        mImageOut=InfoImage(1).Width;
        nImageOut=InfoImage(1).Height;

        %Load image into matlab variable "Image"
        image3D=zeros(nImageOut,mImageOut,NumberImages,'uint16'); %
        for i=1:NumberImages*Number_zStacks
           image3D(:,:,i)=imread(strcat(parent_directory,'\',listing(i).name));
        end
        image3DOut = image3D;
    end
   
[turquoiseStackAll,mImageTurq,nImageTurq] = searchandload(parent_directory,NumberImages,Number_zStacks,'Turquoise');
[venusStackAll,mImageVenus,nImageVenus] = searchandload(parent_directory,NumberImages,Number_zStacks,'Venus');
[draq5StackAll,mImagedraq5,nImagedraq5] = searchandload(parent_directory,NumberImages,Number_zStacks,'DRAQ5');   


turquoiseZStack = zeros(nImageTurq,mImageTurq,NumberImages); %declare
venusZStack = zeros(nImageVenus,mImageVenus,NumberImages); %declare
draq5ZStack = zeros(nImagedraq5,mImagedraq5,NumberImages); %declare


Waitbar_Handle1 = waitbar(0,'Please wait...','Name','SCIENCE in PROGRESS...');
for k = 1:NumberImages:(Number_zStacks*NumberImages)
        
    N=1;
    if(k < Number_zStacks*NumberImages)
        for j = k:(k+NumberImages-1)
            
            turquoiseZStack(:,:,N) = turquoiseStackAll(:,:,j);
            venusZStack(:,:,N) = venusStackAll(:,:,j);
            draq5ZStack(:,:,N) = draq5StackAll(:,:,j);
            
            N = N+1;
        end
    end
    
    fret_rescaling_reslicing_colormap_functions( turquoiseZStack, venusZStack, draq5ZStack, NumberImages, scale, frequencyX, frequencyY, frequencyZ, resliced_directory, smoothing, iterations, steps, reduc, threshold)

end


close all
disp('Done!')
close(Waitbar_Handle1);
        beep on
        beep
        pause(0.5)
        beep
        beep off
end