C++ ArrayXd::size方法代码示例

    inline ArrayXd lm::Dplus(const ArrayXd& d) {
	ArrayXd   di(d.size());
	double  comp(d.maxCoeff() * threshold());
	for (int j = 0; j < d.size(); ++j) di[j] = (d[j] < comp) ? 0. : 1./d[j];
	m_r          = (di != 0.).count();
	return di;
    }

ArrayXd GoSUM::CModelVariables::hcPoint2ModelPoint(const ArrayXd &x)
{
    if ( x.size()!=mvs.size() )   throw "GoSUM::CModelVariables::hcPoint2ModelPoint error: wrong dimension";
    int j,dim=int(x.size());
    ArrayXd X(dim);

    for ( j=0; j<dim; j++ )
    {   X(j)=mvs[j].generateSampleValue(x(j));  }

    return X;
}

void CMATLAB::matPut(string filename,const ArrayXd &X,string Xname)
{
    MATFile *pmat=matOpen(filename,string("w"));
    if (!pmat) throw "CMATLAB::exportTo error: matOpen failed";
    mxArray *pa=mxCreateDoubleMatrix((int)X.size(),1);
    if (!pa) throw "CMATLAB::exportTo error: mxCreateDoubleMatrix failed";
    memcpy((void *)(mxGetPr(pa)), (void *)X.data(), X.size()*sizeof(double));
    if (!matPutVariable(pmat,string("X"),pa)) throw "CMATLAB::exportTo error: matlab.matPutVariable failed";
    mxDestroyArray(pa);
    if (!matClose(pmat)) throw "CMATLAB::exportTo error: matlab.matClose failed";

}

SSFPSimple::SSFPSimple(const ArrayXd &flip, const double TR, const ArrayXd &phi) :
    SteadyState()
{
    m_TR = TR;
    m_flip = (flip * M_PI / 180.).replicate(phi.rows(), 1);
    m_nphi = phi.size();
    m_phi = ArrayXd::Zero(m_flip.size());
    int start = 0;
    for (int i = 0; i < phi.size(); i++) {
        m_phi.segment(start, flip.size()).setConstant(phi[i] * M_PI / 180.);
        start += flip.size();
    }
}

NOMAD::Point CMADS::ArrayXd2NOMADPoint(const ArrayXd &x)
{
    int i,n=int(x.size());
    NOMAD::Point p(n);
    for ( i=0; i<n; i++ ) p[i]=x(i);
    return p;
}

void GoSUM::CModelVariables::setNTuple(const ArrayXd &X,int _at)
{
    if ( X.size()!=mvs.size() ) throw "GoSUM::CModelVariables::setNTuple error: bad nTupe size";

    int i,N=int(mvs.size());
    for ( i=0; i<N; i++ ) mvs[i].setSampleValue(X(i),_at);
}

void CMT::HistogramNonlinearity::setParameters(const ArrayXd& parameters) {
	if(parameters.size() != mHistogram.size())
		throw Exception("Wrong number of parameters.");

	for(int i = 0; i < mHistogram.size(); ++i)
		mHistogram[i] = parameters[i];
}

    const ArrayXd glmDist::devResid(const ArrayXd &y, const ArrayXd &mu, const ArrayXd &wt) const {
	int n = mu.size();
	return as<ArrayXd>(::Rf_eval(::Rf_lang4(as<SEXP>(d_devRes),
						as<SEXP>(NumericVector(y.data(), y.data() + n)),
						as<SEXP>(NumericVector(mu.data(), mu.data() + n)),
						as<SEXP>(NumericVector(wt.data(), wt.data() + n))
					 ), d_rho));
    }

    //@{
    double                   gammaDist::aic     (const ArrayXd& y, const ArrayXd& n, const ArrayXd& mu,
						 const ArrayXd& wt, double dev) const {
	double   nn(wt.sum());
	double disp(dev/nn);
	double   ans(0), invdisp(1./disp);
	for (int i = 0; i < mu.size(); ++i)
	    ans += wt[i] * ::Rf_dgamma(y[i], invdisp, mu[i] * disp, true);
	return -2. * ans + 2.;
    }

double computeBinWidth(const MatrixXd& positions) {
	// assumes first col of positions corresponds to dominant eigenvect
	ArrayXd firstCol = positions.col(0).array();
	firstCol -= firstCol.mean();
	double SSE = firstCol.matrix().squaredNorm();
	double variance = SSE / firstCol.size();
	double std = sqrt(variance);
	double targetBinsPerStd = (MAX_HASH_VALUE - HASH_VALUE_OFFSET) / TARGET_HASH_SPREAD_STDS;
	return std / targetBinsPerStd;
}

    //@{
    double                binomialDist::aic     (const ArrayXd& y, const ArrayXd& n, const ArrayXd& mu,
						 const ArrayXd& wt, double dev) const {
	ArrayXd    m((n > 1).any() ? n : wt);
	ArrayXd   yy((m * y).unaryExpr(Round<double>()));
	m = m.unaryExpr(Round<double>());
	double ans(0.);
	for (int i=0; i < mu.size(); ++i)
	    ans += (m[i] <= 0. ? 0. : wt[i]/m[i]) * ::Rf_dbinom(yy[i], m[i], mu[i], true);
	return (-2. * ans);
    }

ArrayXd GoSUM::CModelVariables::expandNTuple(const ArrayXd &X) const
{
    if ( X.size()!=mvs.size() )  throw "GoSUM::CModelVariables::expand error: wrong X size";
    int i,j,k,N=int(mvs.size()),eN=expandedSize(),exsize;
    ArrayXd eX=ArrayXd::Zero(eN);
    for ( i=j=0; i<N; i++,j+=exsize )
    {   exsize=mvs[i].expandedSize();
        if ( exsize==1 )  {  eX(j)=X(i);  }
        else              {  for ( k=0; k<exsize; k++ )  if (k==X(i)) { eX(j+k) = 1.; break; } } }
    return eX;
}

    void merPredD::updateXwts(const ArrayXd& sqrtXwt) {
        if (d_Xwts.size() != sqrtXwt.size())
            throw invalid_argument("updateXwts: dimension mismatch");
        std::copy(sqrtXwt.data(), sqrtXwt.data() + sqrtXwt.size(), d_Xwts.data());
        if (sqrtXwt.size() == d_V.rows()) { // W is diagonal
            d_V              = d_Xwts.asDiagonal() * d_X;
            for (int j = 0; j < d_N; ++j)
                for (MSpMatrixd::InnerIterator Utj(d_Ut, j), Ztj(d_Zt, j);
                     Utj && Ztj; ++Utj, ++Ztj)
                    Utj.valueRef() = Ztj.value() * d_Xwts.data()[j];
        } else {
            SpMatrixd      W(d_V.rows(), sqrtXwt.size());
            const double *pt = sqrtXwt.data();
            W.reserve(sqrtXwt.size());
            for (Index j = 0; j < W.cols(); ++j, ++pt) {
                W.startVec(j);
                W.insertBack(j % d_V.rows(), j) = *pt;
            }
            W.finalize();
            d_V              = W * d_X;
            SpMatrixd      Ut(d_Zt * W.adjoint());
            if (Ut.cols() != d_Ut.cols())
                throw std::runtime_error("Size mismatch in updateXwts");

            // More complex code to handle the pruning of zeros
            MVec(d_Ut.valuePtr(), d_Ut.nonZeros()).setZero();
            for (int j = 0; j < d_Ut.outerSize(); ++j) {
                MSpMatrixd::InnerIterator lhsIt(d_Ut, j);
                for (SpMatrixd::InnerIterator  rhsIt(Ut, j); rhsIt; ++rhsIt, ++lhsIt) {
                    Index                         k(rhsIt.index());
                    while (lhsIt && lhsIt.index() != k) ++lhsIt;
                    if (lhsIt.index() != k)
                        throw std::runtime_error("Pattern mismatch in updateXwts");
                    lhsIt.valueRef() = rhsIt.value();
                }
            }
        }
        d_VtV.setZero().selfadjointView<Eigen::Upper>().rankUpdate(d_V.adjoint());
        updateL();
    }

    double glmDist::aic(const ArrayXd& y, const ArrayXd& n, const ArrayXd& mu,
			const ArrayXd& wt, double dev) const {
	int nn = mu.size();
	double ans =
	    ::Rf_asReal(::Rf_eval(::Rf_lang6(as<SEXP>(d_aic),
					     as<SEXP>(NumericVector(y.data(), y.data() + nn)),
					     as<SEXP>(NumericVector(n.data(), n.data() + nn)),
					     as<SEXP>(NumericVector(mu.data(), mu.data() + nn)),
					     as<SEXP>(NumericVector(wt.data(), wt.data() + nn)),
					     PROTECT(::Rf_ScalarReal(dev))), d_rho));
	UNPROTECT(1);
	return ans;
    }

    int gesdd(MatrixXd& A, ArrayXd& S, MatrixXd& Vt) {
	int info, mone = -1, m = A.rows(), n = A.cols();
	std::vector<int> iwork(8 * n);
	double wrk;
	if (m < n || S.size() != n || Vt.rows() != n || Vt.cols() != n)
	    throw std::invalid_argument("dimension mismatch in gesvd");
	F77_CALL(dgesdd)("O", &m, &n, A.data(), &m, S.data(), A.data(),
			 &m, Vt.data(), &n, &wrk, &mone, &iwork[0], &info);
	int lwork(wrk);
	std::vector<double> work(lwork);
	F77_CALL(dgesdd)("O", &m, &n, A.data(), &m, S.data(), A.data(),
			 &m, Vt.data(), &n, &work[0], &lwork, &iwork[0], &info);
	return info;
    }