%Program for visualizing SIGMA month files data (Kaupo Komsaare, kaupoko@gmail.com)

clear all
close all

% ----- Those values must be change ------------------------

datafile = 'S1A110200.xl' % SIGMA monthfile
day1 = 2;     % first day of month
day2 = 10;     % last day of month

% -----------------------------------------------------------

%%%%%%%%%%%%%%%%%%%%%%%%%%% Fraction diameter borders %%%%%%%

dpb = [0.422 0.562 0.75 1.00 1.33 1.78 2.37 3.16 4.22 5.62 7.50]';
msize = size (dpb);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

dy = 1; % Number of days plotting on one graph (NOTE: change Xtick properties line 61 to avoid overlapping)

fid=fopen(datafile,'r');
tra=fgets(fid);
tra=fgets(fid);
tra=fgets(fid);
i=0;
disp (' PLEASE WAIT... loading data... ');
while feof(fid)==0
datas=fgets(fid); 
i=i+1;
rawline(i,:)=str2num(datas(1:length(datas)));%

year_(i)=[round(rawline(i,1)/10000)+2000];
month_(i)=[round(mod(rawline(i,1),10000)/100)];
day_(i)=[mod(rawline(i,1),100)];
hour_(i)=[ceil(rawline(i,2)/100)]-1;
minu_(i)=[mod(rawline(i,2),100)];
secu_(i)=0;
time(i)=datenum(year_(i),month_(i),day_(i),hour_(i),minu_(i),secu_(i));

% TE(i)   = rawline(i,4);   % temperature (°C)
% RH(i)   = rawline(i,5);   % relative humidity
% P(i)    = rawline(i,6);   % air pressure (hPa)

CPB(i,:) = [rawline(i,9:18)];    % 10 fractions 
CNB(i,:) = [rawline(i,19:28)];
% ZP(i,:)  = [rawline(i,29:44)];   % 16 fractions
% ZN(i,:)  = [rawline(i,45:60)];
% NP(i)  = [rawline(i,7)];      % noise index +
% NN(i)  = [rawline(i,8)];      % noise index -

end;
fclose(fid);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

disp (' Reading done... ');
[m n] = size(rawline)

CPB(:,11)=NaN; % this empty column is only for plotting purposes
CNB(:,11)=NaN;

year1b = str2num(datafile(4:5)) + 2000;   % year
monthb = str2num(datafile(6:7));        % month

%----------------------- loop for drawing multiple plots ----------------
for day1b=day1:day2 		 
%-------------------------------------------------------------------------

figure
scnsize = get(0,'ScreenSize');
%-----------------------------POSITIVE IONS--------------------------------------------

axes('Position', [0.1 0.55 0.88 0.35]) % plot coordinates [left bottom width height]
pcolor(time,dpb,log10(max(0,CPB'+1e-10))) % draws pcolor plot
set(gca,'FontSize',10,'FontName','Arial','LineWidth',1)% font attributes
shading flat
set(gca,'YScale','log','TickDir','out','XTick',[datenum(year1b,monthb,day1b):2/24:datenum(year1b,monthb,day1b)+dy]) % Y axis is logaritmic, 2/24 means ticks over 2 hours
caxis([0 4]) % scale log10(1)-log10(4)
grid % grid on
axis([datenum(year1b,monthb,day1b)-0.1 datenum(year1b,monthb,day1b)+dy+0.1 0.4 10]) % axis properties (xmin xmax ymin ymax)
datetick('x',15,'keeplimits','keepticks')
ylabel('Size (nm)','FontSize',10) % y label and font size
if dy > 1
    title(['SIGMA(+) Tartu ' datestr(datenum(year1b,monthb,day1b)),' - ',datestr(datenum(year1b,monthb,day1b)+dy-1)],'FontWeight','bold','FontSize',10)
end
if dy == 1
    title(['SIGMA(+) Tartu ' datestr(datenum(year1b,monthb,day1b))],'FontWeight','bold','FontSize',10)
end
set(gca,'YTickLabel',{'1','10'})

% drawing colorbar
cb=colorbar('v6');
set(cb,'Position',[0.9 0.54 0.03 0.38]);
set(cb,'YTickLabel',[1 10 100 1000 10000]);                  
set(cb,'TickLength',[.01 .01]);                        
set(cb,'YLim',[0 4]);                                   
set(cb,'YTick',[0 1 2 3 4]);                             
set(cb,'YColor','k');
set(cb,'FontSize', 10);
set(cb,'FontName','Arial');

%------------------------------NEGATIVE IONS--------------------------------------------
axes('Position', [0.1 0.09 0.88 0.35])
pcolor(time,dpb,log10(max(0,CNB'+1e-10)))
set(gca,'FontSize',10,'FontName','Arial','LineWidth',1)
shading flat
set(gca,'YScale','log','TickDir','out','XTick',[datenum(year1b,monthb,day1b):2/24:datenum(year1b,monthb,day1b)+dy])
caxis([0 4])
grid
axis([datenum(year1b,monthb,day1b)-0.1 datenum(year1b,monthb,day1b)+dy+0.1 0.4 10])
datetick('x',15,'keeplimits','keepticks') %15 is time and 7 is day number in X-axis
xlabel('Time of day','FontSize',12)
ylabel('Size (nm)','FontSize',10)
if dy > 1
    title(['SIGMA(-) Tartu ' datestr(datenum(year1b,monthb,day1b)),' - ',datestr(datenum(year1b,monthb,day1b)+dy-1)],'FontWeight','bold','FontSize',10)
end
if dy == 1
    title(['SIGMA(-) Tartu ' datestr(datenum(year1b,monthb,day1b))],'FontWeight','bold','FontSize',10)
end
set(gca,'YTickLabel',{'1','10'})

% drawing colorbar
cb=colorbar('v6');
set(cb,'Position',[0.9 0.08 0.03 0.38]);
set(cb,'YTickLabel',[1 10 100 1000 10000]);                  
set(cb,'TickLength',[.01 .01]);                        
set(cb,'YLim',[0 4]);                                   
set(cb,'YTick',[0 1 2 3 4]);                             
set(cb,'YColor','k');
set(cb,'FontSize', 10);
set(cb,'FontName','Arial');

% adding text
h = axes('Position',[0 0 1 1],'Visible','off');
set(gcf,'CurrentAxes',h)
text(0.87,0.95,'dN/d(logD)','Rotation',0,'FontSize',10,'FontName','Times New Roman')

%----------------------------- SAVING PLOT into current folder ----------
eval(['print -dpng plot_SIGMA_' datestr(datenum(year1b,monthb,day1b),10) datestr(datenum(year1b,monthb,day1b),5) datestr(datenum(year1b,monthb,day1b),7) '.png'])
close all;

end %   for day1b