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


%ipart = 5 ;    %  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 = ['SZSim2OutRun' num2str(irun)] ;




%vidat = 2 ;
idat = 2 ;
%vin = 3 ;
in = 3 ;
%vchoosesimflag = [1 2 4] ;
simflag = 1 ;


figure(1) ;
clf ;


%  Set up colorful color map
cocomap = [0,    0,   1; ...
          .35, .35, .35; ...
          .5,    0,  .5; ...
           1,    0,   0; ...
           1,   .5,   0; ...
         .35,  .35, .35; ...
           0,    1,   0; ...
           0,    1,   1] ; 
colormap(cocomap) ;



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 ;



n = vn(in) ;


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


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

arraymap = cellMap{simflag+1} ;
     %  3-d array of SiZer maps
     %    simflag + 1, since simflag is 0,1,2,3,4
arrayNRB = cellNRB{simflag+1} ;
     %  2-d matrix of row-wse numbers of pixels Red or Blue


vNRB = sum(arrayNRB,2) ;
    %  Vector of whole image numbers of pixels Red or Blue


[temp,vind] = sort(vNRB) ;
    %  gives pointers to indices in original data 
    %      (ie to case numbers)


%  Make "representative" SiZer Maps
%
%vi = 400:50:1000 ;
vi = [500 750 850 950] ;
    %  indices of plot to make
iplot = 1 ;
for i = vi ;

  subplot(4,1,iplot) ;

  mapout = arraymap(:,:,vind(i)) ;
      %  sizermap at this quantile

  image([left,right],[log10(hmin),log10(hmax)],mapout') ;
    set(gca,'YDir','normal') ;
    if iplot == 1 ;
      title('Figure 2a') ;
      set(gca,'XtickLabel',[]) ;
    elseif iplot == 2 ;
      title('Figure 2b') ;
      set(gca,'XtickLabel',[]) ;
    elseif iplot == 3 ;
      title('Figure 2c') ;
      set(gca,'XtickLabel',[]) ;
    elseif iplot == 4 ;
      title('Figure 2d') ;
    end ;
%    title([num2str(i) 'th ordered SiZer map, from run ' num2str(irun) ...
%           ', idat = ' num2str(idat) ', n = ' num2str(n), ...
%           ', ' simstr]) ;

  vax = axis ;
  text(vax(1) + 0.1 * (vax(2) - vax(1)), ...
       vax(3) + 0.9 * (vax(4) - vax(3)), ...
       ['# Reds or Blues = ' num2str(vNRB(vind(i)))],'Color','w') ;


  iplot = iplot + 1 ;

end ;    %  of i loop



pstr = 'SiZerDist1Fig2' ;

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




