load('Figure3-5.mat');

% study the signals with MLAT<-5deg
valued_south=find((log10(B_str)>-3.75 | log10(E_str)>-2.75) & Plan_3B3E>0.75 & ...
    log10(B_str)<5 & log10(E_str)<5 & mlat<-5 & Wave_P>90 & Wave_3B3E>90);
wna_bin=90:3:175;
B_power=nan(1,length(wna_bin));
signal_num=nan(1,length(wna_bin));
for i=1:length(wna_bin)
    get=find(Wave_3B3E(valued_south)>wna_bin(i) & Wave_3B3E(valued_south)<wna_bin(i)+3);
    signal_num(i)=length(get);
    B_power(i)=mean(B_str(valued_south(get)))^2*length(get)*0.4;
end

f2=figure(2);
p21=plot(wna_bin+1.5,B_power/sum(B_power),'Linewidth',1)
hold on
p22=plot(wna_bin+1.5,signal_num/sum(signal_num),'Linewidth',1)
xlabel('\theta_k (°)')
ylabel('Normalized distribution')
xlim([90 180])
ylimits=ylim;
line([145 145],ylimits,'Linestyle','--','Color','m','Linewidth',1);
leg1=legend([p21 p22],{'Magnetic energy density','Signal pulse occurrence'},'Location','NorthEast');