disp('Running MATLAB script file VarHeavyTail2figE.m') ;

%    Makes Figure E for the Variable Heavy Tails paper

%    Copied from UNC2001RespSiz1.m
%    with parameters set to:
%        ipart = 42  -  Full Data slopes of Log-Log Complementary CDF
%        vitimeb = 10 ;    %  vector of indices of time blocks

%
%    For making analysis of Response Data, based on UNC traces
%    gathered in April 2001, by Don Smith and Co.
%    and filtered by Felix Hernandez Campos.
%
%    Because of the size of the data, they are summarized as:
%         Summary statistics:
%              rsp01power.xls
%              rsp01morestats.xls  (with trimmed down version in 
%                                           rsp01morestatsraw.xls
%                                   and just data saved in
%                                           rsp01morestatsraw.txt)
%         Quantiles of body of distribution:
%              responsesizequantiles.xls
%         Individual responses:
%              20010420_1300.raw
%                     to
%              20010429_800.raw
%                (21 time blocks)
%
%    For further details, see the files:
%         responsesizetailsnotes.txt
%         daydatekey.txt



itimeb = 10 ;    %  vector of indices of time blocks
                    %  1 - Monday, April 23, 8:00-12:00
                    %  2 - Monday, April 23, 13:00-17:00
                    %  3 - Monday, April 23, 19:30-23:30
                    %  4 - Tuesday, April 24, 8:00-12:00
                    %  5 - Tuesday, April 24, 13:00-17:00
                    %  6 - Tuesday, April 24, 19:30-23:30
                    %  7 - Wednesday, April 25, 8:00-12:00
                    %  8 - Wednesday, April 25, 13:00-17:00
                    %  9 - Wednesday, April 25, 19:30-23:30
                    %  10 - Thursday, April 26, 8:00-12:00
                    %  11 - Thursday, April 26, 13:00-17:00
                    %  12 - Thursday, April 26, 19:30-23:30
                    %  13 - Friday, April 20, 8:00-12:00
                    %  14 - Friday, April 20, 13:00-17:00
                    %  15 - Friday, April 20, 19:30-23:30
                    %  16 - Saturday, April 21, 8:00-12:00
                    %  17 - Saturday, April 21, 13:00-17:00
                    %  18 - Saturday, April 21, 19:30-23:30
                    %  19 - Sunday, April 29, 8:00-12:00
                    %  20 - Sunday, April 29, 13:00-17:00
                    %  21 - Sunday, April 29, 19:30-23:30





figure(1) ;
clf ;


%  Load summary data (gives filenames, etc.)
%
%load('\Research\NetworkData\ResponseSize\DataUNC2001\UNC2001RespSiz1DatSumm.mat') ;
load('..\ResponseSize\DataUNC2001\UNC2001RespSiz1DatSumm.mat') ;


%  Load Quantiles
%
%fid = fopen('\Research\NetworkData\ResponseSize\DataUNC2001\responsesizequantilesRaw.txt','rt') ;
fid = fopen('..\ResponseSize\DataUNC2001\responsesizequantilesRaw.txt','rt') ;
        %  'rt' is for "read only" and "text"
fstr = '%f' ;
[mq,cnt] = fscanf(fid,fstr) ;
fclose('all') ;
mq = reshape(mq,21,1000)' ;
        %  use this for 1000x21 matrix, since reshape works column wise



disp(['    Working on ' daystrcell{itimeb} ', ' timestrcell{itimeb}]) ;



%  first get part that appeared earlier
%  (i.e. for given quantiles)
%
vq = mq(:,virowraw(itimeb)) ;
    %  unpack vector of quantiles

n = vnumresp(itimeb) ;
    %  unpack full sample size

nq = length(vq) ;
vccdf = flipud((1/n) + (0:(nq-1))' / nq) ;


%  load large data values, and integrate
%
filename = filenamecell{itimeb} ;
%mtailresp = load(['\Research\NetworkData\ResponseSize\DataUNC2001\' filename]) ;
mtailresp = load(['..\ResponseSize\DataUNC2001\' filename]) ;
vtailsize = mtailresp(:,1) ;
ntail = length(vtailsize) ;

vq = [vq; sort(vtailsize)] ;
vccdf = [vccdf; flipud((1:ntail)' / n)] ;

[vq, vindsort] = sort(vq) ;
vccdf = vccdf(vindsort) ;
nq = length(vq) ;

flagsame = vq(1:(nq-1)) == vq(2:nq) ;
    %  one where next entry is the same
vq = vq(~flagsame) ;
vccdf = vccdf(~flagsame) ;
    %  keep only entries that are different (always keeping last one)
nq = length(vq) ;


%  compute direct difference quotients
%
sbottom = -4.01 ;
stop = 0.01 ;

vcent = mean([log10(vq(2:nq)), log10(vq(1:(nq-1)))],2) ;
vdenom = (log10(vq(2:nq)) - log10(vq(1:(nq-1)))) ;
flag0 = vdenom <= 0 ;
nflag0 = sum(flag0) ;
if nflag0 > 0 ;
  vdenom(flag0) = ones(nflag0,1) ;
  vslope = (log10(vccdf(2:nq)) - log10(vccdf(1:(nq-1)))) ./ vdenom ;
  vslope(flag0) = sbottom * ones(nflag0,1) ;
else ;
  vslope = (log10(vccdf(2:nq)) - log10(vccdf(1:(nq-1)))) ./ vdenom ;
end ;

left = 1 ;
right = 9 ;
bottom = -7 ;
top = 0.01 ;


%  Make gridded difference quotients
%
del1 = 0.3 ;
del2 = 0.03 ;
    %  grid widths
vgrid1 = (left:del1:right)' ;
vgrid2 = (left:del2:right)' ;
    %  grid points covering full range
minq = log10(min(vq)) ;
maxq = log10(max(vq)) ;
    %  endpoints where have data beyond
vgrid1 = vgrid1(minq < vgrid1) ;
vgrid2 = vgrid2(minq < vgrid2) ;
vgrid1 = vgrid1(vgrid1 < maxq) ;
vgrid2 = vgrid2(vgrid2 < maxq) ;
    %  reduce to grid within the data
vheightgrid1 = interp1(log10(vq),log10(vccdf),vgrid1) ;
vheightgrid2 = interp1(log10(vq),log10(vccdf),vgrid2) ;
    %  interpolate to grid
ngrid1 = length(vheightgrid1) ;
ngrid2 = length(vheightgrid2) ;
vcent1 = mean([vgrid1(2:ngrid1), vgrid1(1:(ngrid1-1))],2) ;
vcent2 = mean([vgrid2(2:ngrid2), vgrid2(1:(ngrid2-1))],2) ;
    %  center points for plotting
vslope1 = (vheightgrid1(2:ngrid1) - vheightgrid1(1:(ngrid1-1))) ./ del1 ;
vslope2 = (vheightgrid2(2:ngrid2) - vheightgrid2(1:(ngrid2-1))) ./ del2 ;
    %  get difference quotients (slopes)



%  Make graphics
%
subplot(2,2,1) ;
  plot(log10(vq),log10(vccdf),'b-','LineWidth',2) ;
    axis([left, right, bottom, top]) ;
%    title(['log-log CCDF, HTTP R. S., n = ' num2str(n)],'FontSize',18) ;
%    xlabel('log_{10}(quantile)','FontSize',18) ;
%    xlabel('log_{10}(quantile)','FontSize',12) ;
    xlabel('log_{10}(x)','FontSize',12) ;
%    ylabel('log_{10}(1 - prob.)','FontSize',18) ;
%    ylabel('log_{10}(1 - prob.)','FontSize',12) ;
    ylabel('log_{10}(P[X > x])','FontSize',12) ;
%      tx = left + 0.1 * (right - left) ;
%      ty = bottom + 0.1 * (top - bottom) ;
%    text(tx,ty,[daystrcell{itimeb} ...
%                ', ' timestrcell{itimeb}],'FontSize',18) ;

%    vprobplot = [0.5, 0.9, 0.999, 0.99999] ;
%    vlq = [] ;
%    for probplot = vprobplot ;
%      ty = log10(1 - probplot) ;
%      [temp,ind2min] = min(abs(vccdf - (1 - probplot))) ;
%      tx = log10(vq(ind2min)) ;
%      vlq = [vlq; tx] ;
%      hold on ;
%        plot([tx; tx],[bottom,top],'g--') ;
%        tx = tx + 0.02 * (right - left) ;
%        text(tx,ty,[num2str(probplot) ' quantile'], ...
%                              'FontSize',15,'Color','g') ;
%      hold off ;
%    end ;
%    nlq = length(vlq) ;

subplot(2,2,3) ;
  plot(vcent,vslope,'r-','LineWidth',2) ;
    axis([left,right, sbottom, stop]) ;
%    title(['Corresponding Slopes of log-log CCDF'],'FontSize',18) ;
%    xlabel('log_{10}(quantile)','FontSize',18) ;
%    xlabel('log_{10}(quantile)','FontSize',12) ;
    xlabel('log_{10}(x)','FontSize',12) ;
%    ylabel('-\alpha','FontSize',18) ;
    ylabel('-\alpha','FontSize',12) ;

    hold on ;
      plot([left; right],[-1; -1],'k--') ;
      plot([left; right],[-2; -2],'k--') ;
%      plot([vlq'; vlq'],[sbottom * ones(1,nlq); stop * ones(1,nlq)],'g--') ;
    hold off ;

subplot(2,2,2) ;
  plot(vcent1,vslope1,'m-','LineWidth',2) ;
    axis([left,right, sbottom, stop]) ;
%    title(['\Delta = ' num2str(del1) ' Slopes of log-log CCDF'],'FontSize',18) ;
%    xlabel('log_{10}(quantile)','FontSize',18) ;
%    xlabel('log_{10}(quantile)','FontSize',12) ;
    xlabel('log_{10}(x)','FontSize',12) ;
%    ylabel('-\alpha','FontSize',18) ;
    ylabel('-\alpha','FontSize',12) ;

    hold on ;
      plot([left; right],[-1; -1],'k--') ;
      plot([left; right],[-2; -2],'k--') ;
%      plot([vlq'; vlq'],[sbottom * ones(1,nlq); stop * ones(1,nlq)],'g--') ;
    hold off ;

subplot(2,2,4) ;
  plot(vcent2,vslope2,'m-','LineWidth',2) ;
    axis([left,right, sbottom, stop]) ;
%    title(['\Delta = ' num2str(del2) ' Slopes of log-log CCDF'],'FontSize',18) ;
%    xlabel('log_{10}(quantile)','FontSize',18) ;
%    xlabel('log_{10}(quantile)','FontSize',12) ;
    xlabel('log_{10}(x)','FontSize',12) ;
%    ylabel('-\alpha','FontSize',18) ;
    ylabel('-\alpha','FontSize',12) ;

    hold on ;
      plot([left; right],[-1; -1],'k--') ;
      plot([left; right],[-2; -2],'k--') ;
%      plot([vlq'; vlq'],[sbottom * ones(1,nlq); stop * ones(1,nlq)],'g--') ;
    hold off ;

%      orient landscape ;
%    print('-dpsc2',['UNC2001RS1allCCDFSfull' num2str(itimeb) '.ps']) ;





pstr = 'VarHeavyTail2figE' ;
  orient portrait ;
  set(gcf,'PaperSize',[6.0, 4.5]) ; 
  set(gcf,'PaperPosition',[0.25, 0.25, 5.5, 4.0]) ; 
print('-dpsc', [pstr '.ps']) ;




