C++ MatrixXf::data方法代码示例

/**
 * Create the augmented knots vector and fill in matrix Xt with the spline basis vectors
 */
void basis_splines_endreps_local_v2(float *knots, int n_knots, int order, int *boundaries, int n_boundaries,  Eigen::MatrixXf &Xt) {
  assert(n_boundaries == 2 && boundaries[0] < boundaries[1]);
  int n_rows = boundaries[1] - boundaries[0]; // this is frames so evaluate at each frame we'll need
  int n_cols = n_knots + order;  // number of basis vectors we'll have at end 
  Xt.resize(n_cols, n_rows); // swapped to transpose later. This order lets us use it as a scratch space
  Xt.fill(0);
  int n_std_knots = n_knots + 2 * order;
  float std_knots[n_std_knots];

  // Repeat the boundary knots on there to ensure linear out side of knots
  for (int i = 0; i < order; i++) {
    std_knots[i] = boundaries[0];
  }
  // Our original boundary knots here
  for (int i = 0; i < n_knots; i++) {
    std_knots[i+order] = knots[i];
  }
  // Repeat the boundary knots on there to ensure linear out side of knots
  for (int i = 0; i < order; i++) {
    std_knots[i+order+n_knots] = boundaries[1];
  }
  // Evaluate our basis splines at each frame.
  for (int i = boundaries[0]; i < boundaries[1]; i++) {
    int idx = -1;
    // find index such that i >= knots[idx] && i < knots[idx+1]
    float *val = std::upper_bound(std_knots, std_knots + n_std_knots - 1, 1.0f * i);
    idx = val - std_knots - 1;
    assert(idx >= 0);
    float *f = Xt.data() + i * n_cols + idx - (order - 1); //column offset
    basis_spline_xi_v2(std_knots, n_std_knots, order, idx, i, f);
  }
  // Put in our conventional format where each column is a basis vector
  Xt.transposeInPlace();
}

void VertexBufferObject::update(const Eigen::MatrixXf& M)
{
  assert(id != 0);
  glBindBuffer(GL_ARRAY_BUFFER, id);
  glBufferData(GL_ARRAY_BUFFER, sizeof(float)*M.size(), M.data(), GL_DYNAMIC_DRAW);
  rows = M.rows();
  cols = M.cols();
  check_gl_error();
}

	void PlotBorderPixels(Eigen::MatrixXf& mat, const Graph& graph, WeightMap weights, BorderPixelMap border_pixels)
	{
		float* p = mat.data();
		for(auto eid : as_range(boost::edges(graph))) {
			float v = weights[eid];
			for(unsigned int pid : boost::get(border_pixels, eid)) {
				p[pid] = v;
			}
		}
	}

void run ()
{

  const double tol = argument[2];

#ifdef MRTRIX_UPDATED_API

  const Eigen::MatrixXf in1 = load_matrix<float> (argument[0]);
  const Eigen::MatrixXf in2 = load_matrix<float> (argument[1]);
  
  if (in1.rows() != in2.rows() || in1.cols() != in2.cols())
    throw Exception ("matrices \"" + Path::basename (argument[0]) + "\" and \"" + Path::basename (argument[1]) + "\" do not have matching sizes"
                     " (" + str(in1.rows()) + "x" + str(in1.cols()) + " vs " + str(in2.rows()) + "x" + str(in2.cols()) + ")");
  
  const size_t numel = in1.size();
  for (size_t i = 0; i != numel; ++i) {
    if (std::abs (*(in1.data()+i) - *(in2.data()+i)) > tol)
      throw Exception ("matrices \"" + Path::basename (argument[0]) + "\" and \"" + Path::basename (argument[1]) + " do not match within specified precision of " + str(tol)
                       + " (" + str(*(in1.data()+i)) + " vs " + str(*(in2.data()+i)) + ")");
  }
  
#else

  Math::Matrix<float> in1, in2;
  in1.load (argument[0]);
  in2.load (argument[1]);
  
  if (in1.rows() != in2.rows() || in1.columns() != in2.columns())
    throw Exception ("matrices \"" + Path::basename (argument[0]) + "\" and \"" + Path::basename (argument[1]) + "\" do not have matching sizes"
                     " (" + str(in1.rows()) + "x" + str(in1.columns()) + " vs " + str(in2.rows()) + "x" + str(in2.columns()) + ")");
            
  const size_t numel = in1.rows() * in1.columns();
  for (size_t i = 0; i != numel; ++i) {
    if (std::abs (*(in1.ptr()+i) - *(in2.ptr()+i)) > tol)
      throw Exception ("matrices \"" + Path::basename (argument[0]) + "\" and \"" + Path::basename (argument[1]) + " do not match within specified precision of " + str(tol)
                       + " (" + str(*(in1.ptr()+i)) + " vs " + str(*(in2.ptr()+i)) + ")");
  }               

#endif

  CONSOLE ("data checked OK");
}

void load( Archive & ar,
           Eigen::MatrixXf & t,
           const unsigned int file_version )
{
    int n0;
    ar >> BOOST_SERIALIZATION_NVP( n0 );
    int n1;
    ar >> BOOST_SERIALIZATION_NVP( n1 );
    t.resize( n0, n1 );
    ar >> make_array( t.data(), t.rows()*t.cols() );
}

void save( Archive & ar,
           const Eigen::MatrixXf & t,
           const unsigned int file_version )
{
    int n0 = t.rows();
    ar << BOOST_SERIALIZATION_NVP( n0 );
    int n1 = t.cols();
    ar << BOOST_SERIALIZATION_NVP( n1 );
    ar << boost::serialization::make_array( t.data(),
                                            t.rows()*t.cols() );
}

void run ()
{

  const Eigen::MatrixXf in1 = load_matrix<float> (argument[0]);
  const Eigen::MatrixXf in2 = load_matrix<float> (argument[1]);
  
  if (in1.rows() != in2.rows() || in1.cols() != in2.cols())
    throw Exception ("matrices \"" + Path::basename (argument[0]) + "\" and \"" + Path::basename (argument[1]) + "\" do not have matching sizes"
                     " (" + str(in1.rows()) + "x" + str(in1.cols()) + " vs " + str(in2.rows()) + "x" + str(in2.cols()) + ")");

  const size_t numel = in1.size();

  auto opt = get_options ("frac");
  if (opt.size()) {
  
    const double tol = opt[0][0];
    
    for (size_t i = 0; i != numel; ++i) {
      if (std::abs ((*(in1.data()+i) - *(in2.data()+i)) / (0.5 * (*(in1.data()+i) + *(in2.data()+i)))) > tol)
        throw Exception ("matrices \"" + Path::basename (argument[0]) + "\" and \"" + Path::basename (argument[1]) + " do not match within fractional precision of " + str(tol)
                         + " (" + str(*(in1.data()+i)) + " vs " + str(*(in2.data()+i)) + ")");
    }
  
  } else {
  
    double tol = 0.0;
    opt = get_options ("abs");
    if (opt.size())
      tol = opt[0][0];

    for (size_t i = 0; i != numel; ++i) {
      if (std::abs (*(in1.data()+i) - *(in2.data()+i)) > tol)
        throw Exception ("matrices \"" + Path::basename (argument[0]) + "\" and \"" + Path::basename (argument[1]) + " do not match within absolute precision of " + str(tol)
                         + " (" + str(*(in1.data()+i)) + " vs " + str(*(in2.data()+i)) + ")");
    }
        
  }
  
  CONSOLE ("data checked OK");
}

void FiffStreamThread::sendRawBuffer(Eigen::MatrixXf m_matRawData)
{
    if(m_bIsSendingRawBuffer)
    {
//        qDebug() << "Send RawBuffer to client";

        m_qMutex.lock();

        FiffStream t_FiffStreamOut(&m_qSendBlock, QIODevice::WriteOnly);
        t_FiffStreamOut.write_float(FIFF_DATA_BUFFER,m_matRawData.data(),m_matRawData.rows()*m_matRawData.cols());

        m_qMutex.unlock();

    }
//    else
//    {
//        qDebug() << "Send RawBuffer is not activated";
//    }
}

void
readLDRFile(const std::string& file, Eigen::MatrixXf& data)
{
    std::ifstream in(file.c_str(), std::ios::in | std::ios::binary);
    in.seekg(0, std::ios::end);
    int size = in.tellg();
    in.seekg(0, std::ios::beg);
    
    int num_floats = size / (sizeof(float) / sizeof (char));
    int num_rows = num_floats / 6;
    data.resize(6, num_rows);

    float* row_arr = new float[num_floats];
    in.read((char*)(row_arr), size);

    float* data_arr = data.data();
    for (int k = 0; k < num_floats; k++)
        data_arr[k] = row_arr[k];

    data.transposeInPlace();

    in.close();
}

int main(int argc, char* argv[])
{
    try
    {
        // Select a device and display arrayfire info
        int device = argc > 1 ? std::atoi(argv[1]) : 0;
        af::setDevice(device);
        af::info();

        Eigen::MatrixXf C = Eigen::MatrixXf::Random(1e4, 50); // host array
        af::array in(1e4, 50, C.data());                      // copy host data to device

        af::array u, s_vec, vt;
        svd(u, s_vec, vt, in);
    }
    catch (af::exception& e)
    {
        std::cout << e.what() << '\n';

        return 1;
    }

    return 0;
}

void CloudAnalyzer2D::examinePointEvidence() {
  pointEvidence.clear();

  for (int r = 0; r < R->size(); ++r) {
    Eigen::MatrixXf total = Eigen::MatrixXf::Zero(newRows, newCols);

#pragma omp parallel for schedule(static)
    for (int i = 0; i < total.cols(); ++i) {
      for (int j = 0; j < total.rows(); ++j) {
        for (int k = 0; k < numZ; ++k) {
          const Eigen::Vector3d coord(i, j, k);
          const Eigen::Vector3i src =
              (R->at(r) * (coord - newZZ) + zeroZero)
                  .unaryExpr([](auto v) { return std::round(v); })
                  .cast<int>();

          if (src[0] < 0 || src[0] >= numX || src[1] < 0 || src[1] >= numY ||
              src[2] < 0 || src[2] >= numZ)
            continue;

          if (pointInVoxel->at(src))
            ++total(j, i);
        }
      }
    }

    double average, sigma;
    const float *dataPtr = total.data();
    std::tie(average, sigma) = place::aveAndStdev(
        dataPtr, dataPtr + total.size(), [](auto v) { return v; },
        [](auto v) -> bool { return v; });

    cv::Mat heatMap(newRows, newCols, CV_8UC1, cv::Scalar::all(255));
    for (int j = 0; j < heatMap.rows; ++j) {
      uchar *dst = heatMap.ptr<uchar>(j);
      for (int i = 0; i < heatMap.cols; ++i) {
        const double count = total(j, i);
        if (count > 0) {
          const int gray = cv::saturate_cast<uchar>(
              255.0 *
              ((count - average - sigma) / (3.0 * sigma) - 0.0 * sigma));
          dst[i] = 255 - gray;
        }
      }
    }

    if (FLAGS_preview && doors->size()) {
      std::cout << "Number of doors: " << doors->size() << std::endl;
      cv::Mat out;
      cv::cvtColor(heatMap, out, CV_GRAY2BGR);
      cv::Mat_<cv::Vec3b> _out = out;
      for (auto &d : *doors) {
        Eigen::Vector3d bl = d.corner * FLAGS_scale;
        bl[1] *= -1.0;
        bl = R->at(r).inverse() * bl;

        Eigen::Vector3d xAxis = d.xAxis;
        xAxis[1] *= -1.0;
        xAxis = R->at(r).inverse() * xAxis;

        Eigen::Vector3d zAxis = d.zAxis;
        zAxis[1] *= -1.0;
        zAxis = R->at(r).inverse() * zAxis;

        double w = d.w * FLAGS_scale;

        auto color = randomColor();
        for (int i = 0; i < w; ++i) {
          Eigen::Vector3d pix =
              (bl + i * xAxis + i * zAxis + newZZ).unaryExpr([](auto v) {
                return std::round(v);
              });

          if (pix[0] < 0 || pix[0] >= out.cols || pix[1] < 0 ||
              pix[1] >= out.rows)
            continue;

          _out(pix[1], pix[0]) = color;
        }
      }
      cv::rectshow(out);
    }

    pointEvidence.push_back(heatMap);

    if (FLAGS_preview)
      cv::rectshow(heatMap);
  }
}

	//input vectors as columns (Eigen defaults to column major storage)
	void calcMeanAndCovarWeighed(const Eigen::MatrixXf &input, const Eigen::VectorXd &inputWeights, Eigen::MatrixXf &out_covMat, Eigen::VectorXf &out_mean)
	{
		calcMeanAndCovarWeighed(input.data(),inputWeights.data(),input.rows(),input.cols(),out_covMat,out_mean);

	}

	//input vectors as columns (Eigen defaults to column major storage)
	void calcMeanWeighed(const Eigen::MatrixXf &input, const Eigen::VectorXd &inputWeights, Eigen::VectorXf &out_mean)
	{
		AALTO_ASSERT1(!input.IsRowMajor);
		calcMeanWeighed(input.data(),inputWeights.data(),input.rows(),input.cols(),out_mean);
	}

	// Compute d/df a^T*K*b
	Eigen::MatrixXf kernelGradient( const Eigen::MatrixXf & a, const Eigen::MatrixXf & b ) const {
		Eigen::MatrixXf g = 0*f_;
		lattice_.gradient( g.data(), a.data(), b.data(), a.rows() );
		return g;
	}

	Eigen::MatrixXf PointDensityImpl(const std::vector<T>& seeds, const Eigen::MatrixXf& target, Fx fx, Fy fy)
	{
		// radius of box in which to average cluster density
		constexpr int RHO_R = 3;
		constexpr float cMagicSoftener = 0.5f; // 0.62f;
		constexpr float cRangeLossMult = 1.001f;
		constexpr float cRange = 1.3886f;
		const int rows = target.rows();
		const int cols = target.cols();
		// range of kernel s.t. 99.9% of mass is covered
		Eigen::MatrixXf density = Eigen::MatrixXf::Zero(rows, cols);
		for(const T& s : seeds) {
			const int sx = std::round(fx(s));
			const int sy = std::round(fy(s));
			// compute point density as average over a box
			float rho_sum = 0.0f;
			unsigned int rho_num = 0;
			for(int i=-RHO_R; i<=+RHO_R; ++i) {
				for(int j=-RHO_R; j<=+RHO_R; ++j) {
					const int sxj = sx + j;
					const int syi = sy + i;
					if( 0 <= sxj && sxj < rows &&
						0 <= syi && syi < cols)
					{
						rho_sum += target(sxj, syi);
						rho_num ++;
					}
				}
			}
			if(rho_sum == 0.0f || rho_num == 0) {
				continue;
			}
			const float rho = rho_sum / static_cast<float>(rho_num);
			// seed corresponds to a kernel at position (x,y)
			// with sigma = 1/sqrt(pi*rho)
			// i.e. 1/sigma^2 = pi*rho
			// factor pi is already compensated in kernel
			const float sxf = fx(s);
			const float syf = fy(s);
			const float rho_soft = cMagicSoftener * rho;
			// kernel influence range
			const int R = static_cast<int>(std::ceil(cRange / std::sqrt(rho_soft)));
			const int xmin = std::max<int>(sx - R, 0);
			const int xmax = std::min<int>(sx + R, int(rows) - 1);
			const int ymin = std::max<int>(sy - R, 0);
			const int ymax = std::min<int>(sy + R, int(cols) - 1);
			for(int yi=ymin; yi<=ymax; yi++) {
				for(int xi=xmin; xi<=xmax; xi++) {
					float dx = static_cast<float>(xi) - sxf;
					float dy = static_cast<float>(yi) - syf;
					float d2 = dx*dx + dy*dy;
					float delta = rho_soft * KernelSquare(rho_soft*d2);
					density(xi, yi) += delta;
				}
			}
		}
		{
			const float* psrc = target.data();
			const float* psrc_end = psrc + rows*cols;
			float* pdst = density.data();
			for(; psrc!=psrc_end; ++psrc, ++pdst) {
				if(*psrc == 0.0f) {
					*pdst = 0.0f;
				}
			}
		}
		return cRangeLossMult * density;
	}