/*
 * MATLAB Compiler: 3.0
 * Date: Mon Sep 20 09:40:14 2004
 * Arguments: "-B" "macro_default" "-O" "all" "-O" "fold_scalar_mxarrays:on"
 * "-O" "fold_non_scalar_mxarrays:on" "-O" "optimize_integer_for_loops:on" "-O"
 * "array_indexing:on" "-O" "optimize_conditionals:on" "-B" "sgl" "-m" "-W"
 * "main" "-L" "C" "-t" "-T" "link:exe" "-h" "libmmfile.mlib" "-W" "mainhg"
 * "libmwsglm.mlib" "recarga_multiple" 
 */
#include "scsii.h"
#include "libmatlbm.h"
#include "libmmfile.h"
static mxArray * _mxarray0_;
static mxArray * _mxarray1_;

static double _array3_[21] = { 0.0, 17.0, 26.0, 33.0, 39.0, 45.0, 50.0,
                               55.0, 60.0, 65.0, 70.0, 75.0, 79.0, 83.0,
                               87.0, 91.0, 94.0, 97.0, 98.0, 99.0, 100.0 };
static mxArray * _mxarray2_;

static double _array5_[21] = { 0.0, 5.0, 10.0, 15.0, 20.0, 25.0, 30.0,
                               35.0, 40.0, 45.0, 50.0, 55.0, 60.0, 65.0,
                               70.0, 75.0, 80.0, 85.0, 90.0, 95.0, 100.0 };
static mxArray * _mxarray4_;

static double _array7_[21] = { 0.0, 2.0, 4.0, 7.0, 9.0, 12.0, 15.0, 19.0,
                               23.0, 27.0, 31.0, 35.0, 40.0, 45.0, 51.0,
                               57.0, 63.0, 70.0, 78.0, 87.0, 100.0 };
static mxArray * _mxarray6_;
static mxArray * _mxarray8_;
static mxArray * _mxarray9_;
static mxArray * _mxarray10_;
static mxArray * _mxarray11_;
static mxArray * _mxarray12_;
static mxArray * _mxarray13_;
static mxArray * _mxarray14_;
static mxArray * _mxarray15_;

void InitializeModule_scsii(void) {
    _mxarray0_ = mclInitializeDouble(0.0);
    _mxarray1_ = mclInitializeDouble(1.0);
    _mxarray2_ = mclInitializeDoubleVector(1, 21, _array3_);
    _mxarray4_ = mclInitializeDoubleVector(1, 21, _array5_);
    _mxarray6_ = mclInitializeDoubleVector(1, 21, _array7_);
    _mxarray8_ = mclInitializeDoubleVector(0, 0, (double *)NULL);
    _mxarray9_ = mclInitializeDouble(481.0);
    _mxarray10_ = mclInitializeDouble(2.54);
    _mxarray11_ = mclInitializeDouble(1000.0);
    _mxarray12_ = mclInitializeDouble(10.0);
    _mxarray13_ = mclInitializeDouble(.2);
    _mxarray14_ = mclInitializeDouble(2.0);
    _mxarray15_ = mclInitializeDouble(.8);
}

void TerminateModule_scsii(void) {
    mxDestroyArray(_mxarray15_);
    mxDestroyArray(_mxarray14_);
    mxDestroyArray(_mxarray13_);
    mxDestroyArray(_mxarray12_);
    mxDestroyArray(_mxarray11_);
    mxDestroyArray(_mxarray10_);
    mxDestroyArray(_mxarray9_);
    mxDestroyArray(_mxarray8_);
    mxDestroyArray(_mxarray6_);
    mxDestroyArray(_mxarray4_);
    mxDestroyArray(_mxarray2_);
    mxDestroyArray(_mxarray1_);
    mxDestroyArray(_mxarray0_);
}

static mxArray * Mscsii(int nargout_,
                        mxArray * timerain,
                        mxArray * CNII,
                        mxArray * evaprate,
                        mxArray * data4);

_mexLocalFunctionTable _local_function_table_scsii
  = { 0, (mexFunctionTableEntry *)NULL };

/*
 * The function "mlfScsii" contains the normal interface for the "scsii"
 * M-function from file "u:\final_recarga6\recarga\scsii.m" (lines 1-50). This
 * function processes any input arguments and passes them to the implementation
 * version of the function, appearing above.
 */
mxArray * mlfScsii(mxArray * timerain,
                   mxArray * CNII,
                   mxArray * evaprate,
                   mxArray * data4) {
    int nargout = 1;
    mxArray * y = NULL;
    mlfEnterNewContext(0, 4, timerain, CNII, evaprate, data4);
    y = Mscsii(nargout, timerain, CNII, evaprate, data4);
    mlfRestorePreviousContext(0, 4, timerain, CNII, evaprate, data4);
    return mlfReturnValue(y);
}

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

/*
 * The function "Mscsii" is the implementation version of the "scsii"
 * M-function from file "u:\final_recarga6\recarga\scsii.m" (lines 1-50). 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 y = scs (timerain, CNII, evaprate, data4)
 */
static mxArray * Mscsii(int nargout_,
                        mxArray * timerain,
                        mxArray * CNII,
                        mxArray * evaprate,
                        mxArray * data4) {
    mexLocalFunctionTable save_local_function_table_
      = mclSetCurrentLocalFunctionTable(&_local_function_table_scsii);
    mxArray * y = NULL;
    mxArray * Smax = NULL;
    mxArray * newtimerain = NULL;
    mxArray * i = NULL;
    mxArray * CNI = NULL;
    mxArray * CNIII = NULL;
    mxArray * ARCI = NULL;
    mxArray * ARCII = NULL;
    mxArray * ARCIII = NULL;
    mxArray * datapts = NULL;
    mxArray * antRain = NULL;
    mxArray * Wo = NULL;
    mxArray * Wi = NULL;
    mxArray * Pi = NULL;
    mclCopyArray(&timerain);
    mclCopyArray(&CNII);
    mclCopyArray(&evaprate);
    mclCopyArray(&data4);
    /*
     * %
     * % Accepts a 2-column array, timerain, with columns for time and ppt
     * % Sd, the storage depth, and evaprate, the rate of evaporation.
     * % It goes through the precipitation table and takes off the storage
     * % depth for the entire array.
     * % The evaporation rate must be in compatible units with Sd and timerain
     * % (e.g., if Sd is in inches and the time step in the timerain table is 
     * % hours, then evaprate must be in inches per hour).
     * %
     * Pi=0; Wi=0; Wo=0;
     */
    mlfAssign(&Pi, _mxarray0_);
    mlfAssign(&Wi, _mxarray0_);
    mlfAssign(&Wo, _mxarray0_);
    /*
     * antRain = 0;
     */
    mlfAssign(&antRain, _mxarray0_);
    /*
     * datapts = size(timerain,1); % Determines the number of data points
     */
    mlfAssign(
      &datapts,
      mlfSize(mclValueVarargout(), mclVa(timerain, "timerain"), _mxarray1_));
    /*
     * ARCIII = [0 17 26 33 39 45 50 55 60 65 70 75 79 83 87 91 94 97 98 99 100];
     */
    mlfAssign(&ARCIII, _mxarray2_);
    /*
     * ARCII = [0 5 10 15 20 25 30 35 40 45 50 55 60 65 70 75 80 85 90 95 100];
     */
    mlfAssign(&ARCII, _mxarray4_);
    /*
     * ARCI = [0 2 4 7 9 12 15 19 23 27 31 35 40 45 51 57 63 70 78 87 100];
     */
    mlfAssign(&ARCI, _mxarray6_);
    /*
     * CNIII = interp1(ARCII,ARCIII,CNII);
     */
    mlfAssign(
      &CNIII,
      mlfInterp1(
        mclVv(ARCII, "ARCII"),
        mclVv(ARCIII, "ARCIII"),
        mclVa(CNII, "CNII"),
        NULL));
    /*
     * CNI=interp1(ARCII,ARCI,CNII);
     */
    mlfAssign(
      &CNI,
      mlfInterp1(
        mclVv(ARCII, "ARCII"), mclVv(ARCI, "ARCI"), mclVa(CNII, "CNII"), NULL));
    /*
     * for i = 1:datapts
     */
    {
        int v_ = mclForIntStart(1);
        int e_ = mclForIntEnd(mclVv(datapts, "datapts"));
        if (v_ > e_) {
            mlfAssign(&i, _mxarray8_);
        } else {
            /*
             * newtimerain(i,1) = timerain(i,1); % time
             * if i < 481
             * antRain = antRain + timerain(i,1);
             * else
             * antRain = antRain + timerain(i,1) - timerain(i-480,1);
             * end
             * %if data4(i,3)==0
             * %if antRain > 5.33
             * %CN = CNIII;
             * %elseif antRain < 3.55
             * %CN = CNI;
             * %else
             * %CN = ((antRain-3.55)/(1.78)*(CNIII-CNII))+CNI;
             * %end
             * %end
             * Smax = 2.54*(1000/CNII-10);
             * if data4(i,3)==1 
             * Wo=Wi;
             * Pi=Pi+timerain(i,2);  %Cumm. precip counter
             * Wi=((Pi-0.2*Smax)^2)/(Pi+0.8*Smax);  %SCS Cumm. Runoff Formula
             * if Pi<0.2*Smax
             * Wi=0;
             * end
             * newtimerain(i,2)=Wi-Wo;  %incr. Runoff
             * end
             * if data4(i,3)==0 
             * Pi=0; Wi=0; Wo=0;
             * newtimerain(i,2) = timerain(i,2);% reset to precip counters to zero
             * end  
             * end
             */
            for (; ; ) {
                mclIntArrayAssign2(
                  &newtimerain,
                  mclIntArrayRef2(mclVa(timerain, "timerain"), v_, 1),
                  v_,
                  1);
                if (mclLtBool(mlfScalar(v_), _mxarray9_)) {
                    mlfAssign(
                      &antRain,
                      mclPlus(
                        mclVv(antRain, "antRain"),
                        mclIntArrayRef2(mclVa(timerain, "timerain"), v_, 1)));
                } else {
                    mlfAssign(
                      &antRain,
                      mclMinus(
                        mclPlus(
                          mclVv(antRain, "antRain"),
                          mclIntArrayRef2(mclVa(timerain, "timerain"), v_, 1)),
                        mclIntArrayRef2(
                          mclVa(timerain, "timerain"), v_ - 480, 1)));
                }
                mlfAssign(
                  &Smax,
                  mclMtimes(
                    _mxarray10_,
                    mclMinus(
                      mclMrdivide(_mxarray11_, mclVa(CNII, "CNII")),
                      _mxarray12_)));
                if (mclEqBool(
                      mclIntArrayRef2(mclVa(data4, "data4"), v_, 3),
                      _mxarray1_)) {
                    mlfAssign(&Wo, mclVv(Wi, "Wi"));
                    mlfAssign(
                      &Pi,
                      mclPlus(
                        mclVv(Pi, "Pi"),
                        mclIntArrayRef2(mclVa(timerain, "timerain"), v_, 2)));
                    mlfAssign(
                      &Wi,
                      mclMrdivide(
                        mclMpower(
                          mclMinus(
                            mclVv(Pi, "Pi"),
                            mclMtimes(_mxarray13_, mclVv(Smax, "Smax"))),
                          _mxarray14_),
                        mclPlus(
                          mclVv(Pi, "Pi"),
                          mclMtimes(_mxarray15_, mclVv(Smax, "Smax")))));
                    if (mclLtBool(
                          mclVv(Pi, "Pi"),
                          mclMtimes(_mxarray13_, mclVv(Smax, "Smax")))) {
                        mlfAssign(&Wi, _mxarray0_);
                    }
                    mclIntArrayAssign2(
                      &newtimerain,
                      mclMinus(mclVv(Wi, "Wi"), mclVv(Wo, "Wo")),
                      v_,
                      2);
                }
                if (mclEqBool(
                      mclIntArrayRef2(mclVa(data4, "data4"), v_, 3),
                      _mxarray0_)) {
                    mlfAssign(&Pi, _mxarray0_);
                    mlfAssign(&Wi, _mxarray0_);
                    mlfAssign(&Wo, _mxarray0_);
                    mclIntArrayAssign2(
                      &newtimerain,
                      mclIntArrayRef2(mclVa(timerain, "timerain"), v_, 2),
                      v_,
                      2);
                }
                if (v_ == e_) {
                    break;
                }
                ++v_;
            }
            mlfAssign(&i, mlfScalar(v_));
        }
    }
    /*
     * y = newtimerain;
     */
    mlfAssign(&y, mclVv(newtimerain, "newtimerain"));
    mclValidateOutput(y, 1, nargout_, "y", "scsii");
    mxDestroyArray(Pi);
    mxDestroyArray(Wi);
    mxDestroyArray(Wo);
    mxDestroyArray(antRain);
    mxDestroyArray(datapts);
    mxDestroyArray(ARCIII);
    mxDestroyArray(ARCII);
    mxDestroyArray(ARCI);
    mxDestroyArray(CNIII);
    mxDestroyArray(CNI);
    mxDestroyArray(i);
    mxDestroyArray(newtimerain);
    mxDestroyArray(Smax);
    mxDestroyArray(data4);
    mxDestroyArray(evaprate);
    mxDestroyArray(CNII);
    mxDestroyArray(timerain);
    mclSetCurrentLocalFunctionTable(save_local_function_table_);
    return y;
}
