absfun.c

ASUFIT－Matlab－全局拟合程序

/*
 * MATLAB Compiler: 2.0.1
 * Date: Tue May 08 21:28:16 2001
 * Arguments: "-B" "sgl" "-m" "-W" "mainhg" "-L" "C" "asufit.m" "absfun.m"
 * "absfunfree.m" "absfunshift.m" "absfunwidth.m" "dispfit.m" "dispfitdisp.m"
 * "dispfitfree.m" "dispfitwidth.m" "seqmodfree.m" "spcfun.m" "spcfun1.m"
 * "atamult.m" "aprecon.m" 
 */
#include "absfun.h"

extern mxArray * absLeftEnd1;
extern mxArray * absLeftEnd2;
extern mxArray * absRightEnd1;
extern mxArray * absRightEnd2;
extern mxArray * amps;
extern mxArray * chisq;
extern mxArray * data;
extern mxArray * fit;
extern mxArray * numTao;
extern mxArray * numTaoFixed;
extern mxArray * pulseWidth;
extern mxArray * shift;
extern mxArray * t;
extern mxArray * taoFixed;
extern mxArray * textBox;
extern mxArray * tzero;

/*
 * The function "Mabsfun" is the implementation version of the "absfun"
 * M-function from file "C:\MATLABR11\asufit24\absfun.m" (lines 1-49). It
 * contains the actual compiled code for that M-function. It is a static
 * function and must only be called from one of the interface functions,
 * appearing below.
 */
/*
 * function f = absfun(params)
 */
static mxArray * Mabsfun(int nargout_, mxArray * params) {
    mxArray * f = mclGetUninitializedArray();
    mxArray * a = mclGetUninitializedArray();
    mxArray * ans = mclInitializeAns();
    mxArray * datatemp = mclGetUninitializedArray();
    mclForLoopIterator iterator_0;
    mxArray * j = mclGetUninitializedArray();
    mxArray * sigma = mclGetUninitializedArray();
    mxArray * space = mclGetUninitializedArray();
    mxArray * tTemp = mclGetUninitializedArray();
    mxArray * taoFit = mclGetUninitializedArray();
    mxArray * z = mclGetUninitializedArray();
    mclValidateInputs("absfun", 1, &params);
    /*
     * 
     * global data
     * global fit
     * global amps
     * global numTao
     * global numTaoFixed
     * global taoFixed
     * global t
     * global pulseWidth
     * global tzero
     * global chisq
     * global shift
     * global textBox
     * global absLeftEnd1
     * global absRightEnd1
     * global absLeftEnd2
     * global absRightEnd2
     * 
     * %This is the function referred to in the call to "leastsq".
     * %It assumes the form:
     * %y(x)=a11*exp(-k1*t)+a12*exp(-k2*t)+. . .+a1n*exp(-kn*t)
     * %     a21*exp(-k1*t)+a22*exp(-k2*t)+. . .+a2n*exp(-kn*t)	
     * %     .
     * %     .
     * %     .
     * %     am1*exp(-k1*t)+am2*exp(-k2*t)+. . .+amn*exp(-kn*t)
     * %and makes available the vector of k values and the matrix
     * %of a values to the script file "fit.m".
     * 
     * tTemp = t(absLeftEnd1:absRightEnd1)-t(absLeftEnd1);
     */
    mlfAssign(
      &tTemp,
      mlfMinus(
        mlfIndexRef(t, "(?)", mlfColon(absLeftEnd1, absRightEnd1, NULL)),
        mlfIndexRef(t, "(?)", absLeftEnd1)));
    /*
     * tTemp=tTemp';
     */
    mlfAssign(&tTemp, mlfCtranspose(tTemp));
    /*
     * taoFit(1:numTaoFixed) = taoFixed.^-1;
     */
    mlfIndexAssign(
      &taoFit,
      "(?)",
      mlfColon(mlfScalar(1.0), numTaoFixed, NULL),
      mlfPower(taoFixed, mlfScalar(-1.0)));
    /*
     * taoFit(numTaoFixed + 1:numTaoFixed + numTao) = params(1:numTao).^-1;
     */
    mlfIndexAssign(
      &taoFit,
      "(?)",
      mlfColon(
        mlfPlus(numTaoFixed, mlfScalar(1.0)),
        mlfPlus(numTaoFixed, numTao),
        NULL),
      mlfPower(
        mlfIndexRef(params, "(?)", mlfColon(mlfScalar(1.0), numTao, NULL)),
        mlfScalar(-1.0)));
    /*
     * shift = tzero-t(absLeftEnd1);
     */
    mlfAssign(&shift, mlfMinus(tzero, mlfIndexRef(t, "(?)", absLeftEnd1)));
    /*
     * sigma = pulseWidth / 2.355;
     */
    mlfAssign(&sigma, mlfMrdivide(pulseWidth, mlfScalar(2.355)));
    /*
     * space=zeros(length(tTemp),numTao+numTaoFixed);
     */
    mlfAssign(
      &space, mlfZeros(mlfLength(tTemp), mlfPlus(numTao, numTaoFixed), NULL));
    /*
     * for j = 1:(numTao + numTaoFixed)
     */
    for (mclForStart(
           &iterator_0, mlfScalar(1.0), mlfPlus(numTao, numTaoFixed), NULL);
         mclForNext(&iterator_0, &j);
         ) {
        /*
         * space(:,j) = exp(-taoFit(j)*tTemp).*(1+erf((tTemp-shift-taoFit(j)*sigma^2)/sqrt(2)/sigma))*exp(taoFit(j)*(shift+taoFit(j)*sigma^2/2));
         */
        mlfIndexAssign(
          &space,
          "(?,?)",
          mlfCreateColonIndex(),
          j,
          mlfMtimes(
            mlfTimes(
              mlfExp(
                mlfMtimes(mlfUminus(mlfIndexRef(taoFit, "(?)", j)), tTemp)),
              mlfPlus(
                mlfScalar(1.0),
                mlfErf(
                  mlfMrdivide(
                    mlfMrdivide(
                      mlfMinus(
                        mlfMinus(tTemp, shift),
                        mlfMtimes(
                          mlfIndexRef(taoFit, "(?)", j),
                          mlfMpower(sigma, mlfScalar(2.0)))),
                      mlfSqrt(mlfScalar(2.0))),
                    sigma)))),
            mlfExp(
              mlfMtimes(
                mlfIndexRef(taoFit, "(?)", j),
                mlfPlus(
                  shift,
                  mlfMrdivide(
                    mlfMtimes(
                      mlfIndexRef(taoFit, "(?)", j),
                      mlfMpower(sigma, mlfScalar(2.0))),
                    mlfScalar(2.0)))))));
    /*
     * end
     */
    }
    /*
     * datatemp = data(absLeftEnd2:absRightEnd2,absLeftEnd1:absRightEnd1)';
     */
    mlfAssign(
      &datatemp,
      mlfCtranspose(
        mlfIndexRef(
          data,
          "(?,?)",
          mlfColon(absLeftEnd2, absRightEnd2, NULL),
          mlfColon(absLeftEnd1, absRightEnd1, NULL))));
    /*
     * a = space\datatemp;
     */
    mlfAssign(&a, mlfMldivide(space, datatemp));
    /*
     * z = space*a;
     */
    mlfAssign(&z, mlfMtimes(space, a));
    /*
     * f = z - datatemp;
     */
    mlfAssign(&f, mlfMinus(z, datatemp));
    /*
     * fit = z';
     */
    mlfAssign(&fit, mlfCtranspose(z));
    /*
     * amps = a';
     */
    mlfAssign(&amps, mlfCtranspose(a));
    /*
     * chisq = trace(f'*f);
     */
    mlfAssign(&chisq, mlfTrace(mlfMtimes(mlfCtranspose(f), f)));
    /*
     * set(textBox,'string',['chisq = ' num2str(chisq)]);
     */
    mclAssignAns(
      &ans,
      mlfNSet(
        0,
        textBox,
        mxCreateString("string"),
        mlfHorzcat(mxCreateString("chisq = "), mlfNum2str(chisq, NULL), NULL),
        NULL));
    /*
     * drawnow;
     */
    mlfDrawnow(NULL);
    mclValidateOutputs("absfun", 1, nargout_, &f);
    mxDestroyArray(a);
    mxDestroyArray(ans);
    mxDestroyArray(datatemp);
    mxDestroyArray(j);
    mxDestroyArray(sigma);
    mxDestroyArray(space);
    mxDestroyArray(tTemp);
    mxDestroyArray(taoFit);
    mxDestroyArray(z);
    return f;
}

/*
 * The function "mlfAbsfun" contains the normal interface for the "absfun"
 * M-function from file "C:\MATLABR11\asufit24\absfun.m" (lines 1-49). This
 * function processes any input arguments and passes them to the implementation
 * version of the function, appearing above.
 */
mxArray * mlfAbsfun(mxArray * params) {
    int nargout = 1;
    mxArray * f = mclGetUninitializedArray();
    mlfEnterNewContext(0, 1, params);
    f = Mabsfun(nargout, params);
    mlfRestorePreviousContext(0, 1, params);
    return mlfReturnValue(f);
}

/*
 * The function "mlxAbsfun" contains the feval interface for the "absfun"
 * M-function from file "C:\MATLABR11\asufit24\absfun.m" (lines 1-49). The
 * feval function calls the implementation version of absfun through this
 * function. This function processes any input arguments and passes them to the
 * implementation version of the function, appearing above.
 */
void mlxAbsfun(int nlhs, mxArray * plhs[], int nrhs, mxArray * prhs[]) {
    mxArray * mprhs[1];
    mxArray * mplhs[1];
    int i;
    if (nlhs > 1) {
        mlfError(
          mxCreateString(
            "Run-time Error: File: absfun Line: 1 Column: "
            "0 The function \"absfun\" was called with mor"
            "e than the declared number of outputs (1)"));
    }
    if (nrhs > 1) {
        mlfError(
          mxCreateString(
            "Run-time Error: File: absfun Line: 1 Column: "
            "0 The function \"absfun\" was called with mor"
            "e than the declared number of inputs (1)"));
    }
    for (i = 0; i < 1; ++i) {
        mplhs[i] = NULL;
    }
    for (i = 0; i < 1 && i < nrhs; ++i) {
        mprhs[i] = prhs[i];
    }
    for (; i < 1; ++i) {
        mprhs[i] = NULL;
    }
    mlfEnterNewContext(0, 1, mprhs[0]);
    mplhs[0] = Mabsfun(nlhs, mprhs[0]);
    mlfRestorePreviousContext(0, 1, mprhs[0]);
    plhs[0] = mplhs[0];
}