% file: get_fil_data.i
%
% "get_" file for general analyses of fil_ produced data sets
%

% - - - - Input Parameters :
% define the data set to be read in
%

% For LRF/RMF work...
if ( 1 > 0) {
  % pha: extracted with fil analysis
  % rmf: new dsd/CXC ones copied into appropriate file names
  % arf: NONE available
  data_set_name = "fil_Narrow33";
  %%data_set_name = "fil_1103";
  %%data_set_name = "fil_CapellaSim";
  data_dir = "/nfs/spectra/d2/GtoHAK/Fil_out";
  data_n_spectra = 4;
  data_nh = 0.0;  % x10^22 for wabs
  zpg = 0.0;
  fil_group = 0;
  data_got_rmf = 0;
  data_got_arf = 0;
  rgbprange=0.1; }
% - - - -

% get rid of previous data sets
% could have been from type 2 PHA...
delete_data([1:4]);
delete_rmf([1:4]);
delete_arf([1:4]);

% point to the MEG/HEG -1/+1 data sets
% the four spectra will be loaded in in this order:
megm=1; megp=2; hegm=3; hegp=4;
%
% don;t set the colors so that models can be plotted in
% different color...
%set_data_color(megm,cmegm);
%set_data_color(megp,cmegp);
%set_data_color(hegm,chegm);
%set_data_color(hegp,chegp);

% get the pha files
()=load_data(data_dir+"/"+data_set_name+"_"+"MEG"+"m1_isis.dat");
()=load_data(data_dir+"/"+data_set_name+"_"+"MEG"+"p1_isis.dat");

% don't trust the errors... especially the zero-count ones.
% so fix 'em.
%%(gdlo, gdhi, gdval, gderr) = get_data_counts (megm);
%%variable gdzs = where(gdval == 0);
%%if (length(gdzs) > 1) {
%%  gderr[gdzs] = 0.1;
%%  put_data_counts (megm, gdlo, gdhi, gdval, gderr);  
%%  }
% Arg!  ISIS won't let me set err to 0.1...
% insists on 1.866025e+00 when value is 0.0  :(

% load/assign the rmf files
if (data_got_rmf == 1) {
  ()=load_rmf(data_dir+"/"+data_set_name+"_"+"MEG"+"m1_rmf.fits");
  ()=load_rmf(data_dir+"/"+data_set_name+"_"+"MEG"+"p1_rmf.fits");
  assign_rmf(1,1);
  assign_rmf(2,2);
}

% Get the arf files (Note that they are named appropriately)
if (data_got_arf == 1) {
  ()=load_arf(data_dir+"/"+data_set_name+"_"+"MEG"+"_-1_garf.fits");
  ()=load_arf(data_dir+"/"+data_set_name+"_"+"MEG"+"_1_garf.fits");
  assign_arf(1,1);
  assign_arf(2,2);
}

if ( fil_group > 1) {
% - - - 
% Group the data!
group_data(1,fil_group);
group_data(2,fil_group);

% OK, fix the errors here...
gds = get_data_counts (megm);
gds.err = sqrt(gds.value);
put_data_counts (megm, gds);  
gds = get_data_counts (megp);
gds.err = sqrt(gds.value);
put_data_counts (megp, gds);  
% - - -
  }

% flux correct it so we can see fluxed spectrum if desired...
if (data_got_arf == 1) {
  flux_corr(1,0);
  flux_corr(2,0);
}

% do it all for the HEG spectra...
if (data_n_spectra > 2){
  ()=load_data(data_dir+"/"+data_set_name+"_"+"HEG"+"m1_isis.dat");
  ()=load_data(data_dir+"/"+data_set_name+"_"+"HEG"+"p1_isis.dat");
  if (data_got_rmf == 1) {
    ()=load_rmf(data_dir+"/"+data_set_name+"_"+"HEG"+"m1_rmf.fits");
    ()=load_rmf(data_dir+"/"+data_set_name+"_"+"HEG"+"p1_rmf.fits");
    assign_rmf(3,3);
    assign_rmf(4,4);
  }
  if (data_got_arf == 1) {
    ()=load_arf(data_dir+"/"+data_set_name+"_"+"HEG"+"_-1_garf.fits");
    ()=load_arf(data_dir+"/"+data_set_name+"_"+"HEG"+"_1_garf.fits");
    assign_arf(3,3);
    assign_arf(4,4); 
  }
  if ( fil_group > 1) {
    group_data(3,fil_group);
    group_data(4,fil_group);
    % OK, fix the errors here...
    gds = get_data_counts (hegm);
    gds.err = sqrt(gds.value);
    put_data_counts (hegm, gds);  
    gds = get_data_counts (hegp);
    gds.err = sqrt(gds.value);
    put_data_counts (hegp, gds);  
  }
  if (data_got_arf == 1) {
    flux_corr(3,0);
    flux_corr(4,0);
  }
}


% Show what we have
list_arf;
list_rmf;
list_data;

use_counts;