function [OUTPUT,vals]=bin_to_contour_rowsplit(pixel_width,pixel_height,width_in_meters,height_in_meters,num_contours,filename,percent_lim,split)
%This MAtlab function was created by:
%  Seth Humphries 10-MAY-2004
%  It if for use with the ATLC program and it imports the binary data that ATLC generates, 
%  and creates a contour plot of this data.
%
%NOTE: for large values of pixel_height and pixel_width this will take a LONG time to process
%
%An example of how to call this program from the matlab command prompt is as follows:
%   [outputmatrix,inputvalues]=bin_to_contour(600,600,0.001,0.001,60,'BLAH.U.bin');
%   where BLAH is the name of whatever image process by ATLC, 0.001 is 1mm, and the data is 600 by 600 pixels
%
%The following is an explanation of the values input by the user:
%pixel_width is the width of the original image in pixels
%pixel_height is the height of the original image in pixels
%width_in_meters is the physical width of the original image in meters
%height_in_meters is the physical height of the original image in meters
%   the units of output in the .bin file from ATLC are in joules/(m*pixel^2) those values
%   must be adjusted for the diference in units (should be joules/(m^3) ) 
%   see conversionH and conversionW (which must be the same for the program to work).
%
%num_countours is the number of contours desired by the user. 
%   they are in even increments between 100% and 1% of maximum power
%
%filename is the file from which to import the data, it must be the .bin output from the ATLC program
%   the extension of filename is altered to .dat and "_Xpercent" (where X is percent_lim) is added to the end of the name.
%   The data is exported as comma seperated values so that the data may be read into other 
%   programs such as excel, sigma plot, etc. At the top of the data file statistic values from the
%   data are placed for the user's information
 
percent=percent_lim/100;
%read in the binary data
fid = fopen(filename,'r');
[INPUT,count] = fread(fid,'double');
fclose(fid);
average=1;
area=1;
 
%find the conversion factor
conversionH=height_in_meters/pixel_height;
conversionW=width_in_meters/pixel_width;
 
%verify that the factor is the same in both directions and if so apply it to the matrix
if (conversionW == conversionH)
    conversionH=1/conversionW;
    INPUT = INPUT .*conversionH;
    average=conversionH;
    area=conversionH;
else
%    error=['The height conversion factor is not the same';'as the width conversion factor. The contour ';'plot is not adjusted for units                          ';]
    error=['The height conversion factor is not the same as the width conversion factor. The contour plot is not adjusted for units                          ']
end;
 
%convert the data from a vector to a matrix of dimensions "pixel_height by pixel_width"
Asize=max(size(INPUT));
maximum= max(INPUT);
one=maximum*percent;
row=1;col=1;q=1;r=1;
sum=0;count=0;p=1;
for (i=1:Asize);    
    if (INPUT(i) > 0.0)        
        sum=sum+INPUT(i);
        count=count+1;
        if (INPUT(i) >= one)
            data(p)=INPUT(i);
            p=p+1;
            if (row >=split)
                top_half(q)=INPUT(i);
                q=q+1;
            else
                bot_half(r)=INPUT(i);
                r=r+1;
            end;
        end;
    end;
    OUTPUT(row,col) = INPUT(i);
    if (col == pixel_width)
        row=mod(row,pixel_height)+1;
    end;
    col=mod(col,pixel_width)+1;
end;
average=average*sum/count;  %average is either 1 or conversionH
area = p/(area*area);  %area is either 1 or conversionH
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% this part of the function calculates how to adjust the 
%                                     axes of the figure to fit the 1percent contour nicely on the plot
% ymin=0;ymax=0;xmin=0;xmax=0;
% i=1;j=1;
% while (ymin == 0 && i<=pixel_width)
%     while (ymin == 0 && j<=pixel_height)
%         if (OUTPUT(i,j) >= one)
%             ymin=i;
%         end;
%         j=j+1;
%     end;
%     j=1;
%     i=i+1;
% end;
% 
% i=1;j=1;
% while (xmin == 0 && j<=pixel_height)
%     while (xmin == 0 && i<=pixel_width)
%         if (OUTPUT(i,j) >= one)
%             xmin=j;
%         end;
%         i=i+1;
%     end;
%     i=1;
%     j=j+1;
% end;
% 
% i=pixel_height;j=1;
% while (ymax == 0 && i>=1)
%     while (ymax == 0 && j<=pixel_height)
%         if (OUTPUT(i,j) >= one)
%             ymax=i;
%         end;
%         j=j+1;
%     end;
%     j=1;
%     i=i-1;
% end;
% 
% j=pixel_width;i=1;
% while (xmax == 0 && j>=1)
%     while (xmax == 0 && i<=pixel_width)
%         if (OUTPUT(i,j) >= one)
%             xmax=j;
%         end;
%         i=i+1;
%     end;
%     i=1;
%     j=j-1;
% end;
 
%caculate where to place the contours 
incr=(1/num_contours);
vals=[percent:incr:1].*maximum;
 
w=0.05*pixel_width;
h=0.05*pixel_height;
%create a new figure and plot a contour on it. It is faster for smaller number of contour lines
fig_id=figure;
contourf(OUTPUT,vals); %vals can be replaced with a number specifying how many contours 
    %to draw, matlab then calculates where to place them.
%set(gca,'YLim',[ymin-h ymax+h]);    
%set(gca,'XLim',[xmin-w xmax+w]);    
fsize=max(size(filename));
colormap(hsv)
xlabel('Distance [cm]');
ylabel('Distance [cm]');
ylim([0 pixel_height]);
xlim([0 pixel_width]);
yticks=0:pixel_height/4:pixel_height;
xticks=0:pixel_width/4:pixel_width;
set(gca,'ytick',yticks);
set(gca,'yticklabel',[-6.0 -3.0 0.0 3.0 6.0]);
set(gca,'xtick',xticks);
set(gca,'xticklabel',[-6.0 -3.0 0.0 3.0 6.0]);
 
% title(filename(1:fsize-4));
tempname=[filename(1:fsize-4),num2str(percent_lim),'percent_contour_rowsplit.fig'];
saveas(fig_id,tempname,'fig'); %save in Matlab's native format ".fig"
tempname=[filename(1:fsize-4),num2str(percent_lim),'percent_contour_rowsplit.jpg'];
saveas(fig_id,tempname,'jpg'); %save as a Jpeg image
 
 
% Output the data to a data file as described below
% format will be comma seperate values. this takes SOME time
tempname=filename(fsize-4:fsize);
if (tempname(1:4) ~= '.dat')
    outputname=[filename(1:fsize-4),num2str(percent_lim),'percent_rowsplit.dat'];
    surf_file=[filename(1:fsize-4),num2str(percent_lim),'percent_rowsplit_for_plotting_in_surfer.dat'];
else
    outputname=[filename(1:fsize-4),'unbinary',num2str(percent_lim),'percent_rowsplit.dat'];
    surf_file=[filename(1:fsize-4),'unbinary',num2str(percent_lim),'percent_rowsplit_for_plotting_in_surfer.dat'];
end;
outID=fopen(outputname,'w');
surf_ID=fopen(surf_file,'w');
fprintf(outID,'Stats for all the data:\n');
fprintf(outID,'Number of Pixels greater than %g percent = %g',percent_lim,p);
fprintf(outID,'\nArea of Pixels greater than %g percent in meter^2 = %g',percent_lim,area);
fprintf(outID,'\nMaximum= %g',maximum);
fprintf(outID,'\nMinium= %g',min(data));
fprintf(outID,'\nMean= %g',mean(data));
fprintf(outID,'\nMedian= %g',median(data));
fprintf(outID,'\nStandard Deviation= %g',std(data));
fprintf(outID,'\nCovariance= %g',cov(data));
fprintf(outID,'\nCorrelation Coefficients= %g',corrcoef(data));
fprintf(outID,'\n\n'); %2 lines after the information above
 
fprintf(outID,'Stats for top part of the data:\n');
fprintf(outID,'Number of Pixels greater than %g percent = %g',percent_lim,q);
fprintf(outID,'\nArea of Pixels greater than %g percent in meter^2 = %g',percent_lim,q/(conversionH*conversionH));
fprintf(outID,'\nMaximum= %g',max(top_half));
fprintf(outID,'\nMinium= %g',min(top_half));
fprintf(outID,'\nMean= %g',mean(top_half));
fprintf(outID,'\nMedian= %g',median(top_half));
fprintf(outID,'\nStandard Deviation= %g',std(top_half));
fprintf(outID,'\nCovariance= %g',cov(top_half));
fprintf(outID,'\nCorrelation Coefficients= %g',corrcoef(top_half));
fprintf(outID,'\n\n'); %2 lines after the information above
 
fprintf(outID,'Stats for bottom part of the data:\n');
fprintf(outID,'Number of Pixels greater than %g percent = %g',percent_lim,r);
fprintf(outID,'\nArea of Pixels greater than %g percent in meter^2 = %g',percent_lim,r/(conversionH*conversionH));
fprintf(outID,'\nMaximum= %g',max(bot_half));
fprintf(outID,'\nMinium= %g',min(bot_half));
fprintf(outID,'\nMean= %g',mean(bot_half));
fprintf(outID,'\nMedian= %g',median(bot_half));
fprintf(outID,'\nStandard Deviation= %g',std(bot_half));
fprintf(outID,'\nCovariance= %g',cov(bot_half));
fprintf(outID,'\nCorrelation Coefficients= %g',corrcoef(bot_half));
fprintf(outID,'\n\n'); %2 lines after the information above
 
x_per=100*width_in_meters/pixel_width;%cm
y_per=100*height_in_meters/pixel_height;
for (i=1:pixel_height)%rows
    for(j=1:pixel_width)%columns
        if (OUTPUT(i,j) >= one)
            fprintf(outID,'%g,',OUTPUT(i,j)); % use OUTPUT so that the data structure is intact
            fprintf(surf_ID,'%g,%g,%g\n',(j-pixel_width/2)*x_per,(i-pixel_height/2)*y_per,OUTPUT(i,j));
        else
            fprintf(outID,',');
        end;
        fprintf(surf_ID,'%g,%g,%g\n',(j-pixel_width/2)*x_per,(i-pixel_height/2)*y_per,0);
    end;
    fprintf(outID,'\n');    %at the end of each row put a carraige return
end;
fclose(outID); %close the output data file
fclose(surf_ID);