vector_tools.hpp

Back to Vectors

pastel/sys/

```#ifndef PASTELSYS_VECTOR_TOOLS_HPP
#define PASTELSYS_VECTOR_TOOLS_HPP

#include "pastel/sys/vector_tools.h"
#include "pastel/sys/vector.h"

#include "pastel/sys/mytypes.h"
#include "pastel/sys/ensure.h"

#include "pastel/sys/math_functions.h"

namespace Pastel
{

template <typename Real, int N, typename Expression>
std::ostream& operator<<(std::ostream& stream,
const VectorExpression<Real, N, Expression>& vector)
{
const integer size = vector.size();

for (integer i = 0;i < size;++i)
{
stream << vector[i] << " ";
}

return stream;
}

template <typename Real, int N>
std::istream& operator>>(std::istream& stream,
Vector<Real, N>& vector)
{
const integer size = vector.size();

for (integer i = 0;i < size;++i)
{
stream >> vector[i];
}

return stream;
}

template <typename Real, int N, typename Expression>
Real sum(const VectorExpression<Real, N, Expression>& x)
{
const integer size = x.size();

Real result(x[0]);

for (integer i = 1;i < size;++i)
{
result += x[i];
}

return result;
}

template <typename Real, int N, typename Expression>
Real product(const VectorExpression<Real, N, Expression>& x)
{
const integer size = x.size();

Real result(x[0]);

for (integer i = 1;i < size;++i)
{
result *= x[i];
}

return result;
}

template <
typename Real,
int N>
class VectorUnitAxis
: public VectorExpression<Real, N, VectorUnitAxis<Real, N> >
{
public:
typedef const VectorUnitAxis StorageType;

VectorUnitAxis(
integer axis,
integer size)
: axis_(axis)
, size_(size)
{
}

Real operator[](integer index) const
{
if (index == axis_)
{
return 1;
}

return 0;
}

integer size() const
{
return size_;
}

bool involves(
const void* memoryBegin,
const void* memoryEnd) const
{
return false;
}

bool evaluateBeforeAssignment(
const void* memoryBegin,
const void* memoryEnd) const
{
return false;
}

private:
const integer axis_;
const integer size_;
};

template <typename Real, int N>
VectorUnitAxis<Real, N> unitAxis(integer index)
{
PASTEL_STATIC_ASSERT(N != Dynamic);
PASTEL_STATIC_ASSERT(N > 0);

PENSURE2(index >= 0 && index < N, index, N);

return Pastel::unitAxis<Real, N>(N, index);
}

template <typename Real, int N>
VectorUnitAxis<Real, N> unitAxis(
integer dimension, integer index)
{
PENSURE_OP(dimension, >, 0);
PENSURE2(N == Dynamic || dimension == N,
dimension, N);
PENSURE2(index >= 0 && index < dimension, index, dimension);

return VectorUnitAxis<Real, N>(index, dimension);
}

template <typename Real, int N, typename Expression>
inline Vector<Real, ModifyN<N, N - 1>::Result> shrink(
const VectorExpression<Real, N, Expression>& that)
{
PASTEL_STATIC_ASSERT(N > 1 || N == Dynamic);

const integer size = that.size();

Vector<Real, ModifyN<N, N - 1>::Result> result(size - 1);
for (int i = 0;i < size - 1;++i)
{
result[i] = that[i];
}

return result;
}

template <typename Real, int N, typename Expression>
inline Vector<Real, ModifyN<N, N - 1>::Result> shrink(
const VectorExpression<Real, N, Expression>& that,
integer index)
{
PASTEL_STATIC_ASSERT(N > 1 || N == Dynamic);

const integer size = that.size();

Vector<Real, ModifyN<N, N - 1>::Result> result(size - 1);
for (integer i = 0;i < index;++i)
{
result[i] = that[i];
}
for  (integer i = index;i < size - 1;++i)
{
result[i] = that[i + 1];
}

return result;
}

template <
typename Real,
int N,
typename Expression>
class VectorExtend
: public VectorExpression<Real, ModifyN<N, N + 1>::Result, VectorExtend<Real, N, Expression> >
{
public:
typedef const VectorExtend& StorageType;

VectorExtend(
const Expression& expression,
integer dataIndex,
const Real& data)
: expression_(expression)
, dataIndex_(dataIndex)
, data_(data)
{
PENSURE2(dataIndex >= 0 && dataIndex <= expression.size(),
dataIndex, expression.size());
}

Real operator[](integer index) const
{
if (index < dataIndex_)
{
return expression_[index];
}

if (dataIndex_ < index)
{
return expression_[index - 1];
}

return data_;
}

integer size() const
{
return expression_.size() + 1;
}

bool involves(
const void* memoryBegin,
const void* memoryEnd) const
{
return expression_.involves(memoryBegin, memoryEnd);
}

bool evaluateBeforeAssignment(
const void* memoryBegin,
const void* memoryEnd) const
{
// This is a non-trivial expression.
return expression_.involves(memoryBegin, memoryEnd);
}

private:
typename Expression::StorageType expression_;
const integer dataIndex_;
const Real data_;
};

template <typename Real, int N, typename Expression>
inline VectorExtend<Real, N, Expression> extend(
const PASTEL_NO_DEDUCTION(Real)& left,
const VectorExpression<Real, N, Expression>& right)
{
return VectorExtend<Real, N, Expression>(
(const Expression&)right, 0, left);
}

template <typename Real, int N, typename Expression>
inline VectorExtend<Real, N, Expression> extend(
const VectorExpression<Real, N, Expression>& left,
const PASTEL_NO_DEDUCTION(Real)& right)
{
return VectorExtend<Real, N, Expression>(
(const Expression&)left, left.size(), right);
}

template <typename Real, int N, typename Expression>
inline VectorExtend<Real, N, Expression> extend(
const VectorExpression<Real, N, Expression>& left,
const PASTEL_NO_DEDUCTION(Real)& right,
integer index)
{
return VectorExtend<Real, N, Expression>(
(const Expression&)left, index, right);
}

template <typename Real, int N,
typename Expression>
inline Real dot(
const VectorExpression<Real, N, Expression>& that)
{
return sum(that * that);
}

template <typename Real, int LeftN, int RightN,
typename LeftExpression, typename RightExpression>
inline Real dot(
const VectorExpression<Real, LeftN, LeftExpression>& left,
const VectorExpression<Real, RightN, RightExpression>& right)
{
PASTEL_STATIC_ASSERT(LeftN == Dynamic || RightN == Dynamic ||
LeftN == RightN);

return sum(left * right);
}

template <typename Real, int N, typename Expression,
typename NormBijection>
Real norm2(const VectorExpression<Real, N, Expression>& that,
const NormBijection& normBijection)
{
Real result = 0;

const integer n = that.dimension();
for (integer i = 0;i < n;++i)
{
result,
normBijection.signedAxis(that[i]));
}

return result;
}

template <typename Real, int N, typename Expression>
PASTEL_ENABLE_IF_C(N > 1 || N == Dynamic, Real)
norm(const VectorExpression<Real, N, Expression>& that)
{
return std::sqrt(dot(that, that));
}

template <typename Real, int N, typename Expression>
PASTEL_ENABLE_IF_C(N == 1, Real)
norm(const VectorExpression<Real, N, Expression>& that)
{
return mabs(that[0]);
}

template <typename Real, int N, typename Expression>
Real manhattanNorm(const VectorExpression<Real, N, Expression>& that)
{
return sum(mabs(that));
}

template <typename Real, int N, typename Expression>
Real powerSum(
const VectorExpression<Real, N, Expression>& that,
const PASTEL_NO_DEDUCTION(Real)& metric)
{
PENSURE_OP(metric, >=, 1);
return sum(pow(mabs(that), metric));
}

template <typename Real, int N, typename Expression>
Real minkowskiNorm(
const VectorExpression<Real, N, Expression>& that,
const PASTEL_NO_DEDUCTION(Real)& metric)
{
PENSURE_OP(metric, >=, 1);

return pow(powerSum(that, metric), inverse(metric));
}

template <typename Real, int N, typename Expression>
Real maxNorm(const VectorExpression<Real, N, Expression>& that)
{
return max(mabs(that));
}

template <typename Real, int N>
Vector<Real, N> normalize(const Vector<Real, N>& that)
{
return that / norm(that);
}

template <typename Real, int N, typename Expression>
Vector<Real, N> normalize(
const VectorExpression<Real, N, Expression>& that)
{
return Pastel::normalize(evaluate(that));
}

template <typename Real, int N, typename Expression>
Vector<Real, 2> cross(
const VectorExpression<Real, N, Expression>& that)
{
PASTEL_STATIC_ASSERT(N == 2 || N == Dynamic);
PENSURE_OP(that.size(), ==, 2);

Vector<Real, 2> result(-that[1], that[0]);

return result;
}

template <typename Real, int N, typename ExpressionX,
typename ExpressionY>
Vector<Real, 3> cross(
const VectorExpression<Real, N, ExpressionX>& x,
const VectorExpression<Real, N, ExpressionY>& y)
{
PASTEL_STATIC_ASSERT(N == 3 || N == Dynamic);
PENSURE_OP(x.size(), ==, 3);
PENSURE_OP(y.size(), ==, 3);

Vector<Real, 3> result(
x[1] * y[2] - x[2] * y[1],
x[2] * y[0] - x[0] * y[2],
x[0] * y[1] - x[1] * y[0]);

return result;
}

}

#endif
```