Re: fits in channel space

On Aug 3, 2009, at 2:25 PM, Daniel Dewey wrote:

> Maurice,
>
>  There is an ASCII input file format that load_data will
> read (!).  Do
>
> isis>  help load_data

As Dan says, that will do the read.  ISIS doesn't really care about  
the units so much.  It's just vectors of numbers as far as it's  
concerned.  *However*, any fitting functions you use will care about  
the units.  If you want to use the XSPEC fit functions and pretend  
channel # = keV, you'll have to make sure that you pretend the input  
units are keV as well and make sure that ISIS reads them that way  
(which I think means making sure that the ASCII keyword is there for  
the unit of the X-axis).

However, the default for ISIS if not told otherwise is to pretend that  
things are in Angstrom.  In which case you might want to go ahead and  
treat them that way in plotting and fitting.  Which means that you'll  
need to have some fitting functions in Angstrom.  The ISIS "gauss"  
function is already Angstrom based.  I've also got a few others,  
attached at the end of this, that should get you started.  They are:

% aqpo_fit       : Lorentzian, normalized to QPO rms (counts/bin)
% azfc_fit       : Zero-frequency centered Lorentzian (counts/bin)
% aplaw_fit      : Counts/bin power-law
% abkn_plaw_fit  : Counts/bin broken power-law
% adbkn_plaw_fit : Counts/bin doubly broken power-law
% asinwave_fit   : Counts/bin sine wave
% asqrwave_fit   : Counts/bin square wave

Just drop the "_fit" when using them.  E.g., it's: isis>  
fit_fun("aplaw(1)+aqpo(1)");

  Go ahead an plot them pretending they're Angstroms for the moment.   
It's easy enough to change those axis labels down the road.

-Mike

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%% ANGSTROM COUNTS/BIN FUNCTIONS %%%%%%%%%%%%%%%%% 
%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%
%% QPO (Lorentzian) %%
%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%

% Define a s-lang qpo that works better. Properly averages over bin,
% and thus yields a smoother fit.

define aqpo_fit(lo,hi,par)
{
    variable l,r,q,f,al,ah;

    al = hi;
    ah = lo;

    q = par[1];
    f = par[2];

    % In this formulation, par[0] is rms amplitude

    r = par[0]/(0.5 - atan(-2.*q)/PI);

    l = r^2/(al-ah)/PI *
        ( atan(2.*q*(al-f)/f) - atan(2.*q*(ah-f)/f) );

    return l;
}

add_slang_function("aqpo",["norm [rms]","Q [f/FWHM]","f [Angstrom]"]);

define aqpo_defaults(i)
{
    switch(i)
    {case 0:
       return (0.1,0,0.,1.e8);
    }
    {case 1:
       return (1,0,0.01,1.e3);
    }
    {case 2:
       return (1,0,0.,1.e4);
    }
}

set_param_default_hook("aqpo","aqpo_defaults");

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Zero Freq. Centered Lorentzian  %%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

define azfc_fit(lo,hi,par)
{
    variable a,f,al,ah,l;

    % Go from Angstrom to keV (which we pretend is Fourier Hz),
    % Hence it's atan(al/f) - atan(ah/f) below...

    al = hi;
    ah = lo;
    f = par[1];

    % In this formulation, par[0] is RMS

    a = par[0]^2 * PI/2./f;

    l = a*f/(al-ah) *  ( atan(al/f) - atan(ah/f) );

    return l;
}

add_slang_function("azfc",["norm [rms]","f [Angstrom]"]);

define azfc_defaults(i)
{
    switch(i)
    {case 0:
       return (0.1,0,0.,1.e8);
    }
    {case 1:
       return (1,0,1.e-6,1.e3);
    }
}

set_param_default_hook("azfc","azfc_defaults");

%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Counts/Bin Powerlaw  %%
%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%

define aplaw_fit(lo,hi,par)
{
    variable a,s,al,ah,l;

    al = hi;
    ah = lo;

    a = par[0];
    s = par[1];

    if(s != -1)
    {
       l = a/(s+1) * (al^(s+1) - ah^(s+1));
    }
    else
    {
       l = a * ( log(al) - log(ah) );
    }

    l = (l/(al - ah));

    return l;
}

add_slang_function("aplaw",["norm","slope"]);

define aplaw_defaults(i)
{
    switch(i)
    {case 0:
       return (0.1,0,0.,1.e8);
    }
    {case 1:
       return (-2,0,-5,5);
    }
}

set_param_default_hook("aplaw","aplaw_defaults");

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Broken Counts/Bin Powerlaw  %%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

define abkn_plaw_fit(lo,hi,par)
{
    variable a,sa,b,sb,al,ah,aa,l;
    variable iwl,iwh,liwl,liwh;

    al = hi;
    ah = lo;
    aa = (al+ah)/2.;
    l = _at_email.domain.hidden
    a = par[0];
    sa = par[1];
    b = par[2];
    sb = par[3];

    iwl = where(aa<= b);
    iwh = where(aa > b);

    liwl = length(iwl);
    liwh = length(iwh);

    if(liwl != 0)
    {
       l[iwl] = a * aa[iwl]^sa;
    }

    if(liwh != 0)
    {
       l[iwh] = a * b^sa * ( aa[iwh]/b )^sb;
    }

    return l;
}

add_slang_function("abkn_plaw",["norm","slope_1","f_break  
[Angstrom]","slope_2"]);

define abkn_plaw_defaults(i)
{
    switch(i)
    {case 0:
       return (0.1,0,0.,1.e8);
    }
    {case 1:
       return (0,0,-5,5);
    }
    {case 2:
       return (1,0,0.,1.e4);
    }
    {case 3:
       return (-2,0,-5,5);
    }
}

set_param_default_hook("abkn_plaw","abkn_plaw_defaults");

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Double Broken Counts/Bin Powerlaw  %%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

define adbkn_plaw_fit(lo,hi,par)
{
    variable a,sa,ba,sb,bb,sc,al,ah,aa,l;
    variable iwl,iwm,iwh,liwl,liwm,liwh;

    al = hi;
    ah = lo;
    aa = (al+ah)/2.;
    l = _at_email.domain.hidden
    a = par[0];
    sa = par[1];
    ba = par[2];
    sb = par[3];
    bb = par[4];
    sc = par[5];

    iwl = where(aa <= ba);
    iwm = where(aa > ba and aa <= bb);
    iwh = where(aa > bb);

    liwl = length(iwl);
    liwm = length(iwm);
    liwh = length(iwh);

    if(liwl != 0)
    {
       l[iwl] = a * aa[iwl]^sa;
    }

    if(liwm != 0)
    {
       l[iwm] = a * ba^sa * ( aa[iwm]/ba )^sb;
    }

    if(liwh != 0)
    {
       l[iwh] = a * ba^sa * (bb/ba)^sb * (aa[iwh]/bb)^sc;
    }

    return l;
}

add_slang_function("adbkn_plaw",["norm","slope_1","f_brk_1 [Angstrom]",
                    "slope_2","f_brk_2 [Angstrom]","slope_3"]);

define adbkn_plaw_defaults(i)
{
    switch(i)
    {case 0:
       return (0.1,0,0.,1.e8);
    }
    {case 1:
       return (0,0,-5,5);
    }
    {case 2:
       return (1,0,0.,1.e4);
    }
    {case 3:
       return (-1,0,-5,5);
    }
    {case 4:
       return (10,0,0.,1.e4);
    }
    {case 5:
       return (-2,0,-5,5);
    }
}

set_param_default_hook("adbkn_plaw","adbkn_plaw_defaults");


%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Counts/bin Sine Wave %%
%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%

define asinwave_fit(lo,hi,par)
{
    variable a,f,phs,al,ah,l;

    al = hi;
    ah = lo;

    a = par[0];
    f = par[1];
    phs = par[2];

    l = a * ( -cos(2*PI*(f*al-phs)) + cos(2*PI*(f*ah-phs)) )/2/PI/f/ 
(al-ah);

    return l;
}

add_slang_function("asinwave",["norm","freq","phase"]);

define asinwave_defaults(i)
{
    switch(i)
    {case 0:
       return (0.1,0,0.,1.e8);
    }
    {case 1:
       return (1,0,0,1.e4);
    }
    {case 2:
       return (0,0,0.,1.);
    }
}

set_param_default_hook("asinwave","asinwave_defaults");

%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Counts/bin Square Wave %%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%

define asqrwave_fit(lo,hi,par)
{
    variable a,c,w,r=_at_email.domain.hidden
    r[*]=0;

    a = par[0];
    c = par[1];
    w = par[2];

    iw = where( lo>= c-w/2. and hi <= c+w/2. );
    if(length(iw) >0)
    {
       r[iw] = a;
    }

    iw = where( lo< c-w/2. and hi <= c+w/2. );
    if(length(iw) >0)
    {
       r[iw] = a*(hi[iw] - c+w/2.)/(hi[iw]-lo[iw]);
    }

    iw = where( lo >= c-w/2. and hi > c+w/2. );
    if(length(iw) >0)
    {
       r[iw] = a*(c+w/2. - lo[iw])/(hi[iw]-lo[iw]);
    }

    return r;
}

add_slang_function("asqrwave",["norm","center","width"]);

define asqrwave_defaults(i)
{
    switch(i)
    {case 0:
       return (0.1,0,0.,1.e8);
    }
    {case 1:
       return (1,0,0,1.e4);
    }
    {case 2:
       return (0,0,0.,1.);
    }
}

set_param_default_hook("asqrwave","asqrwave_defaults");

----
You received this message because you are
subscribed to the isis-users list.
To unsubscribe, send a message to
isis-users-request_at_email.domain.hiddenwith the first line of the message as:
unsubscribe