NumCpp  2.12.1
A Templatized Header Only C++ Implementation of the Python NumPy Library
trapz.hpp
Go to the documentation of this file.
1
28#pragma once
29
30#include <string>
31
34#include "NumCpp/Core/Shape.hpp"
35#include "NumCpp/Core/Types.hpp"
36#include "NumCpp/NdArray.hpp"
37
38namespace nc
39{
40 //============================================================================
41 // Method Description:
52 template<typename dtype>
53 NdArray<double> trapz(const NdArray<dtype>& inArray, double dx = 1., Axis inAxis = Axis::NONE)
54 {
56
57 const Shape inShape = inArray.shape();
58 switch (inAxis)
59 {
60 case Axis::NONE:
61 {
62 double sum = 0.;
63 for (uint32 i = 0; i < inArray.size() - 1; ++i)
64 {
65 sum += static_cast<double>(inArray[i + 1] - inArray[i]) / 2. + static_cast<double>(inArray[i]);
66 }
67
68 NdArray<double> returnArray = { sum * dx };
69 return returnArray;
70 }
71 case Axis::COL:
72 {
73 NdArray<double> returnArray(inShape.rows, 1);
74 for (uint32 row = 0; row < inShape.rows; ++row)
75 {
76 double sum = 0;
77 for (uint32 col = 0; col < inShape.cols - 1; ++col)
78 {
79 sum += static_cast<double>(inArray(row, col + 1) - inArray(row, col)) / 2. +
80 static_cast<double>(inArray(row, col));
81 }
82
83 returnArray[row] = sum * dx;
84 }
85
86 return returnArray;
87 }
88 case Axis::ROW:
89 {
90 return trapz(inArray.transpose(), dx, Axis::COL);
91 }
92 default:
93 {
94 THROW_INVALID_ARGUMENT_ERROR("Unimplemented axis type.");
95 return {}; // get rid of compiler warning
96 }
97 }
98 }
99
100 //============================================================================
101 // Method Description:
112 template<typename dtype>
113 NdArray<double> trapz(const NdArray<dtype>& inArrayY, const NdArray<dtype>& inArrayX, Axis inAxis = Axis::NONE)
114 {
115 const Shape inShapeY = inArrayY.shape();
116 const Shape inShapeX = inArrayX.shape();
117
118 if (inShapeY != inShapeX)
119 {
120 THROW_INVALID_ARGUMENT_ERROR("input x and y arrays should be the same shape.");
121 }
122
123 switch (inAxis)
124 {
125 case Axis::NONE:
126 {
127 double sum = 0.;
128 for (uint32 i = 0; i < inArrayY.size() - 1; ++i)
129 {
130 const auto dx = static_cast<double>(inArrayX[i + 1] - inArrayX[i]);
131 sum += dx *
132 (static_cast<double>(inArrayY[i + 1] - inArrayY[i]) / 2. + static_cast<double>(inArrayY[i]));
133 }
134
135 NdArray<double> returnArray = { sum };
136 return returnArray;
137 }
138 case Axis::COL:
139 {
140 NdArray<double> returnArray(inShapeY.rows, 1);
141 for (uint32 row = 0; row < inShapeY.rows; ++row)
142 {
143 double sum = 0;
144 for (uint32 col = 0; col < inShapeY.cols - 1; ++col)
145 {
146 const auto dx = static_cast<double>(inArrayX(row, col + 1) - inArrayX(row, col));
147 sum += dx * (static_cast<double>(inArrayY(row, col + 1) - inArrayY(row, col)) / 2. +
148 static_cast<double>(inArrayY(row, col)));
149 }
150
151 returnArray[row] = sum;
152 }
153
154 return returnArray;
155 }
156 case Axis::ROW:
157 {
158 return trapz(inArrayY.transpose(), inArrayX.transpose(), Axis::COL);
159 }
160 default:
161 {
162 THROW_INVALID_ARGUMENT_ERROR("Unimplemented axis type.");
163 return {}; // get rid of compiler warning
164 }
165 }
166 }
167} // namespace nc
#define THROW_INVALID_ARGUMENT_ERROR(msg)
Definition: Error.hpp:37
#define STATIC_ASSERT_ARITHMETIC(dtype)
Definition: StaticAsserts.hpp:39
size_type size() const noexcept
Definition: NdArrayCore.hpp:4524
self_type transpose() const
Definition: NdArrayCore.hpp:4882
const Shape & shape() const noexcept
Definition: NdArrayCore.hpp:4511
A Shape Class for NdArrays.
Definition: Core/Shape.hpp:41
uint32 rows
Definition: Core/Shape.hpp:44
uint32 cols
Definition: Core/Shape.hpp:45
Definition: Cartesian.hpp:40
Axis
Enum To describe an axis.
Definition: Enums.hpp:36
NdArray< dtype > sum(const NdArray< dtype > &inArray, Axis inAxis=Axis::NONE)
Definition: sum.hpp:46
NdArray< double > trapz(const NdArray< dtype > &inArray, double dx=1., Axis inAxis=Axis::NONE)
Definition: trapz.hpp:53
std::uint32_t uint32
Definition: Types.hpp:40