%%% Plots gratings efficiency from GFW
%fname = 'SOLEIL_g1mm_pos2mm_30deg_3_3ps_l100';

function plot_gratings_energy_gfw_SOLEIL_fcn(fname)

phistepdeg=1;
Radius= 25/2; %mm
Length=230; %mm


%% Reads R2 coeffs
%filename=['/oxford/cauchy/' 'R2_' fname '.txt']
filename=['/oxford/cauchy/' fname '_dE' '.gfwout' ]
fid=fopen(filename);
txt=fread(fid);
fclose(fid);

% Component 9 is the sum of components 3-8
% Component 10 is the sliding sum of components 1-8


fnamerep=strrep(fname,'_',' ');
newline=strfind(char(txt'),10);

wavelength_array=zeros(1,180);
efficiencies_array=zeros(180,3,10);
efficiency_phi_weighted=zeros(180,10);
%start_wvl=30;
%end_wvl=170;

inewline=1;
wvl=90;
iwvl=1;
while ((wvl>=0)&&(wvl<180))
    for ipol=1:3
        %           (char(txt(newline(inewline):newline(inewline+1))'))
        vals=str2num((char(txt(newline(inewline):newline(inewline+1))')));
        if (~isempty(vals))
            if (ipol==1)
                wavelength_array(iwvl)=vals(1);
                wvl=vals(1);
                efficiencies_array(iwvl,ipol,:)=vals(2:11);
            else
                efficiencies_array(iwvl,ipol,:)=vals(1:10);
                
            end% if
            inewline=inewline+1;
        else
            wvl=-1;
        end
    end %for
    if (wvl~=0)
        iwvl=iwvl+1;
    end %if
end %while wvl range

mwvl=iwvl-2

figure(301)
clf
for ipol=1:3
    %figure(ipol)
    %clf
    subplot(3,1,ipol)
    hold on
    phi_weights=0;
    for ival=1:8
        if (mod(ival,3)==0)
            LineStyle='-';
        elseif (mod(ival,3)==1)
            LineStyle='--';
        else
            LineStyle=':';
        end
        plot(wavelength_array(1:mwvl),efficiencies_array(1:mwvl,ipol,ival),'Color',[0 1-((ival-1)/7) (ival-1)/7],'LineStyle',LineStyle)
%        efficiency_phi_weighted(1:mwvl,ipol)=efficiency_phi_weighted(1:mwvl,ipol)+efficiencies_array(1:mwvl,ipol,ival)*2*sqrt(Radius^2-(Length*sin((ival-1)*phistepdeg*pi/180))^2);
%        phi_weights=phi_weights+2*sqrt(Radius^2-(Length*sin((ival-1)*phistepdeg*pi/180))^2);
    end %for ival
    efficiency_phi_weighted=efficiency_phi_weighted./phi_weights;
    plot(wavelength_array(1:mwvl),mean(efficiencies_array(1:mwvl,ipol,1:8),3),'Color',[1 0 0],'LineStyle','-')
%    plot(wavelength_array(1:mwvl),efficiency_phi_weighted(1:mwvl,ipol),'Color',[1 0 0],'LineStyle','-')
    
    hold off
    xlabel('theta');
    ylabel('efficiency R^2');
    %    if (ipol==1)
    legend('Phi=0','2','3','4','5','6','7','Phi max','Avg','Location','EastOutside')
    %    end %if
    if (ipol==1)
        title({['Total -- ' fnamerep ]; ['Phi step (deg)= ' num2str(phistepdeg)]})
    elseif (ipol==2)
        title({['P1 -- ' fnamerep ] ; ['Phi step (deg)= ' num2str(phistepdeg)]})
    else
        title({['P2 -- ' fnamerep ] ; ['Phi step (deg)= ' num2str(phistepdeg)]})
    end
end  %for ipol

print('-dpng',['/Users/delerue/travail/svn/sp_upgrade/trunk/data_acquisition/gratings/gratings_efficiency_gfw_' fname '.png']);

%% Reads energy
%filename=['/oxford/cauchy/' 'Energy_' fname '.txt'];
filename=['/oxford/cauchy/'  fname '_Energy' '.gfwout' ]
fid=fopen(filename);
txt=fread(fid);
fclose(fid);

%%% Is the grating period encoded in the filename
gperiodpos=strfind(fname,'gperiod');
mmpos=min(strfind(fname(gperiodpos:end),'mm'));
gperiod=fname(gperiodpos+8:gperiodpos+mmpos-2);
gperiodval=str2num(gperiod);

fnamerep=strrep(fname,'_',' ');
newline=strfind(char(txt'),10);

values_array=zeros(7,300);

inewline=1;
ang=90;
iang=1;
while ((ang>0)&&(ang<180)&&(inewline<length(newline)))
    vals=str2num((char(txt(newline(inewline):newline(inewline+1))')));
    if (~isempty(vals))
        values_array(:,iang)=vals;
        inewline=inewline+1;
    else
        ang=-1;
    end
    iang=iang+1;
end %while wvl range

%Values at 60°
%xxx

mang=iang-2;

figure(302)
clf
axes('position',[.1  .1  .8  .00001 ])
xlim([20 180]);
xlabel('Wavelength [mm]');
set(gca,'YTickLabel','');
set(gca,'YTick',[]);
anglesFreq=[30:30:180];
set(gca,'XTick',anglesFreq);
wvl=gperiodval*(1-cos(anglesFreq*pi/180));
wvl
for iloop=1:length(wvl)
    wvlLabel{iloop}=num2str(wvl(iloop),3);
end
set(gca,'XTickLabel',wvlLabel);

axes('position',[.1  .2  .8  .0001])
xlim([20 180]);
xlabel('Frequency [THz]');
%set(gca,'XDir','reverse')
set(gca,'YTickLabel','');
set(gca,'YTick',[]);
anglesFreq=[30:30:180];
set(gca,'XTick',anglesFreq);
freqs=(3e-1/gperiodval)*(1-cos(anglesFreq*pi/180));
for iloop=1:length(freqs)
    freqsLabel{iloop}=num2str(freqs(iloop),3);
end
set(gca,'XTickLabel',freqsLabel);
%xlim([  (3e-1/gperiodval)*(1-cos(180*pi/180)) ])
xlim([20 180]);
axes('position',[.1  .3  .8  .6])
hold on
grid on
plot(values_array(1,1:mang),values_array(4,1:mang),'-r','LineWidth',3)
plot(values_array(1,1:mang),values_array(6,1:mang),'-.b','LineWidth',2)
plot(values_array(1,1:mang),values_array(7,1:mang),':g','LineWidth',2)
plot(values_array(1,1:mang),values_array(6,1:mang)+values_array(7,1:mang),'o','LineWidth',2)
hold off
xlabel('theta');
ylabel('Energy/bunch (J) (25mm diameter)');
legend('Total','P1','P2','p1+P2','Location','Best');
title({ [fnamerep ] ; 'Expected energy '});
print('-dpng',['/Users/delerue/travail/svn/sp_upgrade/trunk/data_acquisition/gratings/gratings_energy_gfw_' fname '.png']);


end %function

