disp('Running MATLAB script file gpsz1dt.m') ;
%
%    FOR DEVELOPMENT AND TESTING OF MATLAB FUNCTION gpsz1d,
%    General Purpose SiZer for Dependent data,
%    This analyzes time series data


itest = 2 ;     %  1 - sunspot data
                %  2 - lynx data
                %  3 - share yield data
                %  4 - detrended share yield data
                %  5 - Buffalo Snowfall
                %  6 - Australian chocolates data
                %  7 - detrended Australian chocolates data


figure(1) ;
clf ;



if itest == 1 ;    %  Sunspot data

  dtitstr = 'Sunspots' ;
  infstr = 'sunspot.dat' ;
  frmtstr = '%g' ;
  insize = [1, inf] ;
  vicol = [1] ;
  datatran = 'data = sqrt(data) ;' ;
  xstr = 'Year' ;
  ystr = 'sqrt(# Sunspots)' ;

  disp(['    Loading ' dtitstr ' Data']) ;
    infstr = ['\matlab\steve\data\' infstr] ;
    fid = fopen(infstr,'r') ;
      data = fscanf(fid,frmtstr,insize) ;
      data = data' ;
          %  since data rows are read in as columns
      data = data(:,vicol) ;
    fclose(fid) ;
  disp(['    Read in ' num2str(size(data,1)) ' data points']) ;

  if length(datatran) ~= 0 ;    %  Then do some data transformation
    eval(datatran) ;
  end ;


  gpsz1d(data) ;


  subplot(4,4,1) ;
    title('Sunspot Data') ;
    ylabel('sqrt(sunspots)') ;    


    orient landscape ;
  print -dpsc DepSiZerSunspots.ps ;
  


elseif itest == 2 ;    %  Lynx data

  dtitstr = 'Canadian Lynx' ;
           %  Famous time series data, from lynx trappings
  infstr = 'canlynx.dat' ;
  frmtstr = '%g' ;
  insize = [1,inf] ;
  vicol = [1] ;
  datatran = 'data = log10(data + 1) ;' ;
  xstr = 'Year' ;
  ystr = 'log10(# Lynx trapped + 1)' ;

  disp(['    Loading ' dtitstr ' Data']) ;
    infstr = ['\matlab\steve\data\' infstr] ;
    fid = fopen(infstr,'r') ;
      data = fscanf(fid,frmtstr,insize) ;
      data = data' ;
          %  since data rows are read in as columns
      data = data(:,vicol) ;
    fclose(fid) ;
  disp(['    Read in ' num2str(size(data,1)) ' data points']) ;

  if length(datatran) ~= 0 ;    %  Then do some data transformation
    eval(datatran) ;
  end ;


  gpsz1d(data) ;


  subplot(4,4,1) ;
    title('Lynx Trapping Data') ;
    ylabel('log10(lynx + 1)') ;    


    orient landscape ;
  print -dpsc DepSiZerLynx.ps ;



elseif itest == 3 ;    %  share yield data

  dtitstr = 'Share Yield' ;
  %  from Probal Chaudhuri
  %  Index of Stock Profits on "ordinary shares", Indian Stock Market
  %  on scale of "interest rates"
  %  from Reserve Bank of India Bulletin
  %  time series of monthly data, April 1982 - March 1991
  infstr = 'sharyiel.dat' ;
  frmtstr = '%g' ;
  insize = [1, inf] ;
  vicol = [1] ;
  datatran = [] ;
  xstr = 'month' ;
  ystr = 'Share yield' ;

  disp(['    Loading ' dtitstr ' Data']) ;
    infstr = ['\matlab\steve\data\' infstr] ;
    fid = fopen(infstr,'r') ;
      data = fscanf(fid,frmtstr,insize) ;
      data = data' ;
          %  since data rows are read in as columns
      data = data(:,vicol) ;
    fclose(fid) ;
  disp(['    Read in ' num2str(size(data,1)) ' data points']) ;

  gpsz1d(data) ;


  subplot(4,4,1) ;
    title('Share Yield Data') ;


    orient landscape ;
  print -dpsc DepSiZerShareYield.ps ;



elseif itest == 4 ;    %  share yield data

  dtitstr = 'Detrended Share Yield' ;
  %  from Probal Chaudhuri
  %  Index of Stock Profits on "ordinary shares", Indian Stock Market
  %  on scale of "interest rates"
  %  from Reserve Bank of India Bulletin
  %  time series of monthly data, April 1982 - March 1991
  infstr = 'sharyiel.dat' ;
  frmtstr = '%g' ;
  insize = [1, inf] ;
  vicol = [1] ;
  datatran = [] ;
  xstr = 'month' ;
  ystr = 'Share yield' ;

  disp(['    Loading ' dtitstr ' Data']) ;
    infstr = ['\matlab\steve\data\' infstr] ;
    fid = fopen(infstr,'r') ;
      data = fscanf(fid,frmtstr,insize) ;
      data = data' ;
          %  since data rows are read in as columns
      data = data(:,vicol) ;
    fclose(fid) ;
  disp(['    Read in ' num2str(size(data,1)) ' data points']) ;

  n = length(data) ;
  vtime = (1:n)' ;
  vcoeff = polyfit(vtime,data,1) ;
  data = data - polyval(vcoeff,vtime) ;


  gpsz1d(data) ;


  subplot(4,4,1) ;
    title('Detrended Share Yields') ;


    orient landscape ;
  print -dpsc DepSiZerDetShareY.ps ;
  


elseif itest == 5 ;    %  Buffalo Snowfall data


  dtitstr = 'Snowfalls' ;
  infstr = 'buffalo.sno' ;

  disp(['    Loading ' dtitstr ' Data']) ;
    infstr = ['\matlab\steve\data\' infstr] ;
    fid = fopen(infstr,'r') ;
      data = fscanf(fid,frmtstr,insize) ;
      data = data' ;
          %  since data rows are read in as columns
      data = data(:,vicol) ;
    fclose(fid) ;
  disp(['    Read in ' num2str(size(data,1)) ' data points']) ;


  gpsz1d(data) ;


  subplot(4,4,1) ;
    title('Buffalo Snowfalls') ;
    ylabel('inches') ;


    orient landscape ;
  print -dpsc DepSiZerBuffSno.ps ;
  


elseif itest == 6 ;    %  Chocolates data


  load chocs.txt ;
  data = chocs ;


  gpsz1d(data) ;


  subplot(4,4,1) ;
    title('Australian Chocolates') ;


    orient landscape ;
  print -dpsc DepSiZerChocs.ps ;
  


elseif itest == 7 ;    %  deterended chocolates data


  load chocs.txt ;
  data = chocs ;

  n = length(data) ;
  vtime = (1:n)' ;
  vcoeff = polyfit(vtime,data,1) ;
  data = data - polyval(vcoeff,vtime) ;


  gpsz1d(data) ;


  subplot(4,4,1) ;
    title('Detrended Chocolates') ;


    orient landscape ;
  print -dpsc DepSiZerDetChocs.ps ;
  


end ;


