% file: fit_demo.sl    10/28/12 - dd
%
% Do some fitting of the E0102 R+S+B model to the imaging data.
%

% - - - - -
% Setup Data, Model and Viewing parameters

% Load the data and the model
.source e0102_acis_data
;

.source e0102_model
;

% Use coarser cells and blur the data (and model) to be more comparable to the Einstein HRI:
e2d_set_view(1, [ -30 , 30 , 2.5 ] , [ -30 , 30 , 2.5 ] );
e2d_view[1].nsmoo=3;
% set a higher residual color-coding threshold:
e2d_view[1].chigreen=2.5;
% and set a minimum fractional sigma (10%):
e2d_meta[1].scisigma=0.1;


% - - - - -
% The simple scans changed the model parameters to:
s3d_ps.norm1=0.62;
s3d_ps.norm2=0.42;
s3d_ps.norm3=0.028;
%
s3d_ps.az0=-125.0;

% Giving a pre-fitting FOM:
e2d_pars2fom;
s3d_list_model;
e2d_norm_optimum(0.001);
%   - Minimum FOM of 1100.21 at norm factor = 1.0102
e2d_view_resid;


#iffalse
% Comment out the following to manually paste the steps in on
% the command line.

% - - - - -
% Now demonstrate CG fitting of the Ring:
%   xc, yc, Ri, Ro (with Si fixed at Ro), Ra,
%   and the ring norm (different from the overall norm)
% See how fitting is setup by doing:
%   hydra> e2d_cg_fit;
%   ... extensive help info ...

public define e2d_update_pars () 
{ 
   % this line is the same for all cg fits: 
   variable pvals = cg_p0 + cg_smand*cg_dir*e2d_scan_vals[e2d_iscan]; 
   % map the fitting parameter values to the desired model params: 
   s3d_ps.xc = pvals[0];
   s3d_ps.yc = pvals[1];
   s3d_ps.Ri = pvals[2]; 
   s3d_ps.Ro = pvals[3];  s3d_ps.Si = pvals[3];
   s3d_ps.Ra = pvals[4];
   s3d_ps.norm1 = pvals[5];
   % don't allow Ri to be too close to or greater than Ro:
   if ((s3d_ps.Ri+0.6) > s3d_ps.Ro) s3d_ps.Ri = s3d_ps.Ro - 0.6;
};
cg_p0 = [0.0, 3.0, 12.0, 17.0, 56.0, 0.62];
cg_smand = [1.0, 1.0, 1.0, 1.0, 5.0, 0.05];

% do the fit:
e2d_cg_fit(1); 

% deactivate the scanpar definition:
e2d_update_pars_NULL;

% Make some "fomlog" plots:
e2d_fomlog_plot("xc","yc",10.);
e2d_fomlog_plot("Ri","Ro",10.);
e2d_fomlog_plot("Ra","norm1",10.);
e2d_fomlog_plot("Ra","Ro",10.);
e2d_fomlog_history(["xc","yc","Ri","Ro","Ra","norm1"]);


% Set the results:
s3d_ps.xc = -0.30;
s3d_ps.yc = 3.30;
s3d_ps.Ri = 13.2; 
s3d_ps.Ro = 16.60;
  s3d_ps.Si = 16.60;
s3d_ps.Ra = 53.0;
s3d_ps.norm1 = 0.62; 


e2d_pars2fom;
s3d_list_model;
e2d_norm_optimum(0.001);
%   - Minimum FOM of 846.542 at norm factor = 1.03951
e2d_view_resid;

% The FOM has been improved ...

#endif