%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Program to do ratio photobleach correction. % % Need 16-bit individual tiff images in a folder with consecutive numbers % starting with "xxx1" for data up to 9 planes, "xxx01" for up to 99 planes, % and "xxx001" for up to 999 planes. % Plus MS-Excel "XXXX.xls" file with Plane/Numerator/Denominator in 3 % consecutive columns, with first cell in each column containing % descriptions (i.e., "Plane", "I-SO", "GFP" etc. % Plane numbering should start from 1 to however many images. % Intensity data for I-SO plane number 1 should never be Zero. % I-SO intensity can be set to zero in cells after the plane number 1 to % prevent them from being used in the curve fitting, if justified. % % author: Louis Hodgson 2004 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear all; clc; % getting the ratio image file names and intensities file names I=input('Enter number of ratio images: '); % This sets the number of corrections to loop fname=input('Enter file name without running numbers and file extension: ','s'); % This is the name of the ratio file to be corrected dname=input('Enter Excel file name for the intensities, Plane#/Num/Denom; without Excel file extensions: ','s'); % This sets file name for intensities DNAMEIN=sprintf('%s.xls', dname); data_all = xlsread(DNAMEIN); Plane_initial = data_all(1:I); % Plane array goes from 1 to I num_Int = data_all(I+1:2*I); % absolute value intensities for the numerator denom_Int = data_all(2*I+1:3*I); % absolute value intensities for the denominator num_Int_norm = num_Int ./ num_Int(1); % calculates normalized intensity values for numerator, normalized to t=0 denom_Int_norm = denom_Int ./ denom_Int(1); % calculates normalized intensity values for denominator, normalized to t=0 Plane = (Plane_initial - 1).'; % makes plane array start from zero instead of 1 for data fitting purposes Ratio = (num_Int_norm ./ denom_Int_norm).'; % calculates the ratio of the normalized instensities outliers = ~excludedata(Plane, Ratio, 'range', [0 0.1]); % sets fit exclusion range in Y values [fresult,gof] = fit(Plane, Ratio, 'exp2', 'exclude', outliers); % double exponential fit, Plane numbers as X and normalized Ratio as Y display (fresult); % displays fitted parameters and coefficients display (gof); % displays goodness of fit parameters Fit_ratio = fresult(Plane); % producing fitted curve using the plane numbers CF = 1 ./ Fit_ratio; % correction factor is 1 over the ratio decay function a=num2str(fresult.a); b=num2str(fresult.b); c=num2str(fresult.c); d=num2str(fresult.d); rsq=num2str(gof.rsquare); coef=sprintf('a = %s ; b = %s', a, b); coef2=sprintf('c = %s ; d = %s', c, d); coef3=sprintf('r^2 = %s', rsq); plot(Plane, Ratio,'rd', Plane, fresult(0:I-1),'b-'); % Plots results title('Double exponential fit to Photobleach data'); xlabel('Plane Number'); ylabel('Normalized Ratio Decay'); legend('Data', 'Fit'); axis([0, I, 0, 1.2]); text (2, 0.4, 'y = a exp(b*x) + c exp(d*x)'); text (2, 0.33, coef); text (2, 0.26, coef2); text (2, 0.19, coef3); % Writes data files for output and plotting off-line fid=fopen('Fit_Decay_curve.txt', 'a'); fprintf(fid,'%12.6f \r', Fit_ratio); fclose(fid); fid=fopen('Actual_Ratio_data.txt', 'a'); fprintf(fid,'%12.6f \r', Ratio); fclose(fid); fid=fopen('goodness_of_fit.txt', 'a'); fprintf(fid,'%12.6f \r', gof); fclose(fid); % Image processing starting at this point, using the CF as the correction if I < 100; if I < 10; for mult_rep=1:I; counter=mult_rep; FNAMEIN=sprintf('%s%i.tif', fname, counter); ratio=imread(FNAMEIN,'tif'); R=double(ratio); CR=R .* CF(counter); %Array element by element multiplication CR_out=uint16(CR); FNAMEOUT=sprintf('out_%i.tif',counter); imwrite(CR_out, FNAMEOUT, 'tif','compression','none'); % write output file end; else; for mult_rep=1:I; counter=mult_rep; if counter < 10; FNAMEIN=sprintf('%s0%i.tif', fname, counter); else; FNAMEIN=sprintf('%s%i.tif', fname, counter); end; ratio=imread(FNAMEIN,'tif'); R=double(ratio); CR=R .* CF(counter); %Array element by element multiplication CR_out=uint16(CR); if counter < 10; FNAMEOUT=sprintf('out_0%i.tif',counter); else; FNAMEOUT=sprintf('out_%i.tif',counter); end; imwrite(CR_out, FNAMEOUT, 'tif','compression','none'); % write output file end; end; elseif I >= 100; for mult_rep=1:I; counter=mult_rep; if counter < 10; FNAMEIN=sprintf('%s00%i.tif', fname, counter); ratio=imread(FNAMEIN,'tif'); R=double(ratio); CR=R .* CF(counter); %Array element by element multiplication CR_out=uint16(CR); FNAMEOUT=sprintf('out_00%i.tif',counter); imwrite(CR_out, FNAMEOUT, 'tif','compression','none'); % write output file elseif counter < 100; FNAMEIN=sprintf('%s0%i.tif', fname, counter); ratio=imread(FNAMEIN,'tif'); R=double(ratio); CR=R .* CF(counter); %Array element by element multiplication CR_out=uint16(CR); FNAMEOUT=sprintf('out_0%i.tif',counter); imwrite(CR_out, FNAMEOUT, 'tif','compression','none'); % write output file elseif counter >= 100; FNAMEIN=sprintf('%s%i.tif', fname, counter); ratio=imread(FNAMEIN,'tif'); R=double(ratio); CR=R .* CF(counter); %Array element by element multiplication CR_out=uint16(CR); FNAMEOUT=sprintf('out_%i.tif',counter); imwrite(CR_out, FNAMEOUT, 'tif','compression','none'); % write output file end; end; end;