disp('Running MATLAB script file SiZerDist1Fig6.m') ;
%
%    For Making Figure 6
%    in the Advanced Distribution Theory for SiZer paper
%
%
%    Copied from: SiZerSimul2anal.m via SiZerDist1Fig7.m
%
%    Does Analysis of outputs from SiZerSimul2.m
%    which did simulations investigating Simultaneity issues for SiZer
%


%ipart = 2 ;    %  1 - make table of overall percentage of significant pixels
               %  2 - make table of row-wise percetnage of significant pixels
               %  3 - make gray level image of pixel wise percentage (non-gray denom)
               %  4 - make gray level image of pixel wise percentage (nsim denom)
               %  5 - make plots of SiZer maps at specified quantiles
               %             (no links between simul. methods)
               %  6 - make plots of SiZer maps at specified quantiles
               %             (ordered by first simul. method)


irun = 1 ;    %  run number:    1,2
              %    (allows fiddling with parameters)


%
%  Meaning of vsimflag:       %  parameter controlling simulation type
                              %    0  -  Use Pointwise C.I.'s
                              %    1  -  Use old Row-wise Simultaneous C.I.'s
                              %    2 (or not specified) - Use Row-wise based on Hsing at al. (1996)
                              %    3  -  Use global adjustment (sum independent blocks calculated as in 2)
                              %    4  -  Use global adjustment based on Bonferoni and Hsing at al. (1996).
%  Analyze, depending on run
%
if irun == 1 ;

  %  parameters from SiZerSimul2.m
  %  
  vidat = [1 2 3] ;    %  index of data type
                       %  1  -  Uniform(0,1) density estimation
                       %  2  -  N(0,1) regression (equally spaced)
                       %  3  -  N(0,1) regression (Uniform(0,1) design)

  vn = [100 400 1600] ;    %  sample sizes
  
  %  
  %  Parameters loaded from the file:
  %      cellMap cellNRB A2dNNG A2dRM 
  %      vsimflag nsim n left right nbin hmin hmax nh
  %      iregtdist


elseif irun == 2 ;

%  %  parameters from SiZerSimul2.m
%  %
  vidat = [2 3] ;    %  index of data type
                       %  1  -  Uniform(0,1) density estimation
                       %  2  -  N(0,1) regression (equally spaced)
                       %  3  -  N(0,1) regression (Uniform(0,1) design)

  vn = [100 400 1600] ;    %  sample sizes
  
  %  
  %  Parameters loaded from the file:
  %      cellMap cellNRB A2dNNG A2dRM 
  %      vsimflag nsim n left right nbin hmin hmax nh
  %      iregtdist


end ;



infilestr = ['C:\BigDataJSM\SiZerSimulOutputs\SZSim2OutRun' num2str(irun)] ;
%infilestr = ['..\SZSim2OutRun' num2str(irun)] ;


%  Create Data matrices
%
PercSp1Lt1 = [] ;
    %  matrix of percentages, Subplot 1, Linetype 1
PercSp1Lt2 = [] ;
PercSp1Lt3 = [] ;
PercSp2Lt1 = [] ;
PercSp2Lt2 = [] ;
PercSp2Lt3 = [] ;
PercSp3Lt1 = [] ;
PercSp3Lt2 = [] ;
PercSp3Lt3 = [] ;

for idat = vidat ;    %  loop through data types


  if idat == 1 ;
    datstr = 'Uniform(0,1) density estimation' ;
  elseif idat == 2 ;
    datstr = 'N(0,1) regression (equally spaced)' ;
  elseif idat == 3 ;
    datstr = 'N(0,1) regression (Uniform(0,1) design)' ;
  end ;


  for in = 1:length(vn) ;    %  loop through sample sizes

   n = vn(in) ;


    %  load simulation output file
    %
    instr = [infilestr 'Dat' num2str(idat) 'n' num2str(in)] ;
    eval(['load ' instr]) ;
    disp(['    finished read of ' instr]) ;


    lvh = log10(logspace(log10(hmin),log10(hmax),nh)) ;


    %  loop through simultaneity types
    %
    for isimflag = 2:length(vsimflag) ;
               %  Note: starting at 2 skips the pointwise

      simflag = vsimflag(isimflag) ;

      if simflag == 0 ;
        simstr = 'Pointwise C.I.s' ;
      elseif simflag == 1 ;
        simstr = 'Old Row-wise Simultaneous C.I.s' ;
      elseif simflag == 2 ;
        simstr = 'Hsings Row-wise Simultaneous C.I.s' ;
      elseif simflag == 3 ;
        simstr = 'Global Simultaneous C.I.s' ;
      elseif simflag == 4 ;
        simstr = 'Bonferoni Global Simultaneous C.I.s' ;
      end ;


      A2dNRB = cellNRB{isimflag} ;



        flagNRBG = sum(A2dNRB,2) > 0 ;
            %   one for simuls where have a Red or Blue pixel
        percentsig = 100 * sum(flagNRBG) / nsim ;

             %  Global percentage
             %  we do not have to account for non-gray
             %  as all pictures have some non-gray pixels at 
             %  the course scales. 



%        flagNNG = A2dNNG > 0 ;
%            %   one for simuls where have a non-gray pixel 
%            %   as we do not want to count rows that were all gray.
%       vpercentng = 100 * sum(flagNNG,1) / nsim ;
%
%
%        flagNRB = A2dNRB > 0 ;
%            %   one for simuls where have a Red or Blue pixel
%        vpercentsig = 100 * sum(flagNRB,1) ./ (nsim * vpercentng / 100) ;

        if in == 1 ;    %  settings are subplots

          if idat == 1 ;    %  ns are linetypes
            PercSp1Lt1 = [PercSp1Lt1; percentsig] ;
          elseif idat == 2 ;
            PercSp1Lt2 = [PercSp1Lt2; percentsig] ;
          elseif idat == 3 ;
            PercSp1Lt3 = [PercSp1Lt3; percentsig] ;
          end ;

        elseif in == 2 ;    %  settings are subplots

          if idat == 1 ;    %  ns are linetypes
            PercSp2Lt1 = [PercSp2Lt1; percentsig] ;
          elseif idat == 2 ;
            PercSp2Lt2 = [PercSp2Lt2; percentsig] ;
          elseif idat == 3 ;
            PercSp2Lt3 = [PercSp2Lt3; percentsig] ;
          end ;

        elseif in == 3 ;    %  settings are subplots

          if idat == 1 ;    %  ns are linetypes
            PercSp3Lt1 = [PercSp3Lt1; percentsig] ;
          elseif idat == 2 ;
            PercSp3Lt2 = [PercSp3Lt2; percentsig] ;
          elseif idat == 3 ;
            PercSp3Lt3 = [PercSp3Lt3; percentsig] ;
          end ;

        end ;    %  of idat if-block



    end ;    %  of isimflag loop


  end ;    %  of n loop through sample sizes


end ;    %  of idat loop through data types




%  Make main graphic
%
figure(1) ;
clf ;
%left = min(lvh) ;
%right = max(lvh) ;
top = 80 ;

left = 1;
right = 3;
xxx = [1 2 3];

plot(xxx,[PercSp1Lt1(1),PercSp2Lt1(1),PercSp3Lt1(1)],'g:', ...
   xxx,[PercSp1Lt1(2),PercSp2Lt1(2),PercSp3Lt1(2)],'r:', ...
   xxx,[PercSp1Lt1(3),PercSp2Lt1(3),PercSp3Lt1(3)],'b:', ...
   xxx,[PercSp1Lt1(4),PercSp2Lt1(4),PercSp3Lt1(4)],'m:', ...
   xxx,[PercSp1Lt2(1),PercSp2Lt2(1),PercSp3Lt2(1)],'g--', ...
   xxx,[PercSp1Lt2(2),PercSp2Lt2(2),PercSp3Lt2(2)],'r--', ...
   xxx,[PercSp1Lt2(3),PercSp2Lt2(3),PercSp3Lt2(3)],'b--', ...
   xxx,[PercSp1Lt2(4),PercSp2Lt2(4),PercSp3Lt2(4)],'m--', ...
   xxx,[PercSp1Lt3(1),PercSp2Lt3(1),PercSp3Lt3(1)],'g-', ...
   xxx,[PercSp1Lt3(2),PercSp2Lt3(2),PercSp3Lt3(2)],'r-', ...
   xxx,[PercSp1Lt3(3),PercSp2Lt3(3),PercSp3Lt3(3)],'b-', ...
   xxx,[PercSp1Lt3(4),PercSp2Lt3(4),PercSp3Lt3(4)],'m-', ...
   [1; 3],[5; 5],'k-') ;
axis([left,right,0,top]) ;
xlabel('Sample size') ;
ylabel('Global % Sig.') ;
set(gca,'XTick',xxx)
set(gca,'XtickLabel',vn)

title('Figure 6') ;
text(left + 0.05 * (right - left), ...
   0 + 0.9 * (top - 0), ...
   'Conventional','Color','g') ;
text(left + 0.05 * (right - left), ...
   0 + 0.8 * (top - 0), ...
   'Row-Wise','Color','r') ;
text(left + 0.05 * (right - left), ...
   0 + 0.7 * (top - 0), ...
   'Glob. Indep.','Color','b') ;
text(left + 0.05 * (right - left), ...
   0 + 0.6 * (top - 0), ...
   'Glob. Bonf.','Color','m') ;

text(left + 0.7 * (right - left), ...
   0 + 0.85 * (top - 0), ...
   'KDE') ;
text(left + 0.7 * (right - left), ...
   0 + 0.75 * (top - 0), ...
   'FDR') ;
text(left + 0.7 * (right - left), ...
   0 + 0.65 * (top - 0), ...
   'RDR') ;
hold on ;
ljump = 0.0 ;
plot([left + 0.8 * (right - left) ; left + 0.95 * (right - left)], ...
     (0 + (0.85 + ljump) * (top - 0)) * ones(2,1), 'k:') ;
plot([left + 0.8 * (right - left) ; left + 0.95 * (right - left)], ...
     (0 + (0.75 + ljump) * (top - 0)) * ones(2,1), 'k--') ;
plot([left + 0.8 * (right - left) ; left + 0.95 * (right - left)], ...
     (0 + (0.65 + ljump) * (top - 0)) * ones(2,1), 'k-') ;
hold off ;



pstr = 'SiZerDist1Fig6' ;

  orient portrait ;
  set(gcf,'PaperSize',[6.0, 3.0]) ; 
  set(gcf,'PaperPosition',[0.25, 0.25, 5.5, 2.5]) ; 
print('-dpsc', [pstr '.ps']) ;





