meshgrid.c

《精通Matlab6.5版》张志勇等编著

     * %
     * %   For example, to evaluate the function  x*exp(-x^2-y^2) over the 
     * %   range  -2 < x < 2,  -2 < y < 2,
     * %
     * %       [X,Y] = meshgrid(-2:.2:2, -2:.2:2);
     * %       Z = X .* exp(-X.^2 - Y.^2);
     * %       mesh(Z)
     * %
     * %   MESHGRID is like NDGRID except that the order of the first two input
     * %   and output arguments are switched (i.e., [X,Y,Z] = MESHGRID(x,y,z)
     * %   produces the same result as [Y,X,Z] = NDGRID(y,x,z)).  Because of
     * %   this, MESHGRID is better suited to problems in cartesian space,
     * %   while NDGRID is better suited to N-D problems that aren't spatially
     * %   based.  MESHGRID is also limited to 2-D or 3-D.
     * %
     * %   See also SURF, SLICE, NDGRID.
     * 
     * %   J.N. Little 1-30-92, CBM 2-11-92.
     * %   Copyright 1984-2001 The MathWorks, Inc. 
     * %   $Revision: 5.13 $  $Date: 2001/04/15 12:02:40 $
     * 
     * if nargout<3, % 2-D array case
     */
    if (nargout_ < 3) {
        /*
         * if nargin == 1, y = x; end
         */
        if (nargin_ == 1) {
            mlfAssign(&y, mclVsa(x, "x"));
        }
        /*
         * if isempty(x) | isempty(y)
         */
        {
            mxArray * a_ = mclInitialize(mclVe(mlfIsempty(mclVa(x, "x"))));
            if (mlfTobool(a_)
                || mlfTobool(mclOr(a_, mclVe(mlfIsempty(mclVa(y, "y")))))) {
                mxDestroyArray(a_);
                /*
                 * xx = []; yy = []; zz = [];
                 */
                mlfAssign(&xx, _mxarray4_);
                mlfAssign(yy, _mxarray4_);
                mlfAssign(zz, _mxarray4_);
            /*
             * else
             */
            } else {
                mxDestroyArray(a_);
                /*
                 * xx = full(x(:)).'; % Make sure x is a full row vector.
                 */
                mlfAssign(
                  &xx,
                  mlfTranspose(
                    mclVe(
                      mlfFull(
                        mclVe(
                          mclArrayRef1(
                            mclVsa(x, "x"), mlfCreateColonIndex()))))));
                /*
                 * yy = full(y(:));   % Make sure y is a full column vector.
                 */
                mlfAssign(
                  yy,
                  mlfFull(
                    mclVe(
                      mclArrayRef1(mclVsa(y, "y"), mlfCreateColonIndex()))));
                /*
                 * nx = length(xx); ny = length(yy);
                 */
                mlfAssign(&nx, mlfScalar(mclLengthInt(mclVv(xx, "xx"))));
                mlfAssign(&ny, mlfScalar(mclLengthInt(mclVv(*yy, "yy"))));
                /*
                 * xx = xx(ones(ny, 1),:);
                 */
                mlfAssign(
                  &xx,
                  mclArrayRef2(
                    mclVsv(xx, "xx"),
                    mlfOnes(mclVv(ny, "ny"), _mxarray5_, NULL),
                    mlfCreateColonIndex()));
                /*
                 * yy = yy(:,ones(1, nx));
                 */
                mlfAssign(
                  yy,
                  mclArrayRef2(
                    mclVsv(*yy, "yy"),
                    mlfCreateColonIndex(),
                    mlfOnes(_mxarray5_, mclVv(nx, "nx"), NULL)));
            }
        /*
         * end
         */
        }
    /*
     * 
     * else, % 3-D array case
     */
    } else {
        /*
         * if nargin == 1, y = x; z = x; end
         */
        if (nargin_ == 1) {
            mlfAssign(&y, mclVsa(x, "x"));
            mlfAssign(&z, mclVsa(x, "x"));
        }
        /*
         * if nargin ==2, error('Not enough input arguments.'); end
         */
        if (nargin_ == 2) {
            mlfError(_mxarray6_);
        }
        /*
         * if isempty(x) | isempty(y) | isempty(z)
         */
        {
            mxArray * a_ = mclInitialize(mclVe(mlfIsempty(mclVa(x, "x"))));
            if (mlfTobool(a_)) {
                mlfAssign(&a_, mlfScalar(1));
            } else {
                mlfAssign(&a_, mclOr(a_, mclVe(mlfIsempty(mclVa(y, "y")))));
            }
            if (mlfTobool(a_)
                || mlfTobool(mclOr(a_, mclVe(mlfIsempty(mclVa(z, "z")))))) {
                mxDestroyArray(a_);
                /*
                 * xx = []; yy = []; zz = [];
                 */
                mlfAssign(&xx, _mxarray4_);
                mlfAssign(yy, _mxarray4_);
                mlfAssign(zz, _mxarray4_);
            /*
             * else
             */
            } else {
                mxDestroyArray(a_);
                /*
                 * nx = prod(size(x)); ny = prod(size(y)); nz = prod(size(z));
                 */
                mlfAssign(
                  &nx,
                  mlfProd(
                    mclVe(mlfSize(mclValueVarargout(), mclVa(x, "x"), NULL)),
                    NULL));
                mlfAssign(
                  &ny,
                  mlfProd(
                    mclVe(mlfSize(mclValueVarargout(), mclVa(y, "y"), NULL)),
                    NULL));
                mlfAssign(
                  &nz,
                  mlfProd(
                    mclVe(mlfSize(mclValueVarargout(), mclVa(z, "z"), NULL)),
                    NULL));
                /*
                 * xx = reshape(full(x(:)),[1 nx 1]); % Make sure x is a full row vector.
                 */
                mlfAssign(
                  &xx,
                  mlfReshape(
                    mclVe(
                      mlfFull(
                        mclVe(
                          mclArrayRef1(
                            mclVsa(x, "x"), mlfCreateColonIndex())))),
                    mlfHorzcat(_mxarray5_, mclVv(nx, "nx"), _mxarray5_, NULL),
                    NULL));
                /*
                 * yy = reshape(full(y(:)),[ny 1 1]); % Make sure y is a full column vector.
                 */
                mlfAssign(
                  yy,
                  mlfReshape(
                    mclVe(
                      mlfFull(
                        mclVe(
                          mclArrayRef1(
                            mclVsa(y, "y"), mlfCreateColonIndex())))),
                    mlfHorzcat(mclVv(ny, "ny"), _mxarray5_, _mxarray5_, NULL),
                    NULL));
                /*
                 * zz = reshape(full(z(:)),[1 1 nz]); % Make sure z is a full page vector.
                 */
                mlfAssign(
                  zz,
                  mlfReshape(
                    mclVe(
                      mlfFull(
                        mclVe(
                          mclArrayRef1(
                            mclVsa(z, "z"), mlfCreateColonIndex())))),
                    mlfHorzcat(_mxarray5_, _mxarray5_, mclVv(nz, "nz"), NULL),
                    NULL));
                /*
                 * xx = xx(ones(ny,1),:,ones(nz,1));
                 */
                mlfAssign(
                  &xx,
                  mlfIndexRef(
                    mclVsv(xx, "xx"),
                    "(?,?,?)",
                    mlfOnes(mclVv(ny, "ny"), _mxarray5_, NULL),
                    mlfCreateColonIndex(),
                    mlfOnes(mclVv(nz, "nz"), _mxarray5_, NULL)));
                /*
                 * yy = yy(:,ones(1,nx),ones(nz,1));
                 */
                mlfAssign(
                  yy,
                  mlfIndexRef(
                    mclVsv(*yy, "yy"),
                    "(?,?,?)",
                    mlfCreateColonIndex(),
                    mlfOnes(_mxarray5_, mclVv(nx, "nx"), NULL),
                    mlfOnes(mclVv(nz, "nz"), _mxarray5_, NULL)));
                /*
                 * zz = zz(ones(ny,1),ones(nx,1),:);
                 */
                mlfAssign(
                  zz,
                  mlfIndexRef(
                    mclVsv(*zz, "zz"),
                    "(?,?,?)",
                    mlfOnes(mclVv(ny, "ny"), _mxarray5_, NULL),
                    mlfOnes(mclVv(nx, "nx"), _mxarray5_, NULL),
                    mlfCreateColonIndex()));
            }
        /*
         * end
         */
        }
    /*
     * end
     */
    }
    mclValidateOutput(xx, 1, nargout_, "xx", "meshgrid");
    mclValidateOutput(*yy, 2, nargout_, "yy", "meshgrid");
    mclValidateOutput(*zz, 3, nargout_, "zz", "meshgrid");
    mxDestroyArray(nx);
    mxDestroyArray(ny);
    mxDestroyArray(ans);
    mxDestroyArray(nz);
    mxDestroyArray(z);
    mxDestroyArray(y);
    mxDestroyArray(x);
    mclSetCurrentLocalFunctionTable(save_local_function_table_);
    return xx;
}