% Date: 09/23/2003

% Written to run in ISIS (=>v1.1.3), but convertible to a pure s-lang script.
%
% Run time: roughly N^2 routine.  N ~ O(200,000) is limit for an overnight run
%
% Save this and associated files as: 'bblocks_driver.sl', 'sitar.sl', and
% 'sitar_examp.sl".
%
% To run:   isis> .load bblocks_driver  % and then you're finished!
%
% Variables that need to be set.

   variable file="ps.evt.gz";   % Event file to read
   variable ncp_prior = 3.;     % Significance level (see below)
   variable plot_step = 200.;   % Step size (plotting only) of input lightcurve
   variable ctype = 1;          % How to assign intervals to events (see
                                %    sitar_make_data_cells.sl) ctype=3 is for  
                                %    a binned lightcurve
   variable bsize = 1.14104*3.; % binsize for ctype=3. If <=0, then binsize = 
                                %    timedel from FITS file
   variable clump = 0.5;        % Events closer together than this time are 
                                %    assigned to the same cell (see 
                                %    sitar_make_data_cells)
   variable plotit = 0;         % =0 for screen plotting, > 0  print to file
   variable talt = 1;           % start & stop times to use on the lightcurve:
                                %    = 0 for FITS file, > 0 this file (below)
   variable writeit = 0;        % = 0 don't write data, > 0 send to a file

% Load the required subroutines. In ISIS, the first line can (probably should) 
% be replaced by: require("sitar");

   () = evalfile("sitar.sl");
   () = evalfile("bblocks_examp.sl");

%  *** This is the read of the data ***

   variable event_times, tstart, tstop, frame, dtcor, object;
   (event_times, tstart, tstop, frame, dtcor, object) = sitar_examp_read(file);

   variable tstart_file = tstart, tstop_file = tstop;

% Alternative definition for tstart and tstop

   if(talt)
   { 
      tstart = min(event_times);
      tstop = max(event_times);
   }

   if(ctype == 3 and bsize > 0.)
   {
      frame = bsize;
   }

% Make the data cells to be searched.  

   variable cell;
   cell = sitar_make_data_cells( event_times,ctype,clump,frame,tstart,tstop );

% Chandra specific deadtime correction

   cell.dtcor = cell.dtcor * dtcor;

% Main call.  
%
% ncp_prior ~ - natural log of (1 - significance level)
% ncp_prior = 1.1 is ~68% confidence level *per block*
% ncp_prior = 3.0 is ~95%
% ncp_prior = 4.6 is ~99%
% ncp_prior = 7.0 is ~99.9%
%
% ***  This is only a rough way of looking at this parameter which
%      is very bad for ncp_prior low, or for very long lightcurves.  ***

   variable results = sitar_global_optimum( cell, ncp_prior, ctype );

% Bin up the events for plotting purposes

   variable ev;
   ev = sitar_examp_bin( event_times, tstart, tstop, dtcor, plot_step );

% Plot the results (first make sure plotting defaults of linear axes are set,
% and clear any previously defined ranges)

   sitar_examp_plot( 
                    results, ev, object, file, plot_step, ncp_prior,
                    tstart_file, tstop_file, plotit
                   );

% Write the results to a file, if so desired. 

   if(writeit)
   {
      sitar_examp_write( results, object, file, plot_step, ncp_prior );
   }