NumCpp  2.7.0
A Templatized Header Only C++ Implementation of the Python NumPy Library
gradient.hpp
Go to the documentation of this file.
1 #pragma once
29 
34 #include "NumCpp/Core/Shape.hpp"
35 #include "NumCpp/Core/Types.hpp"
36 #include "NumCpp/NdArray.hpp"
37 
38 #include <complex>
39 #include <string>
40 
41 namespace nc
42 {
43  //============================================================================
44  // Method Description:
54  template<typename dtype>
56  {
58 
59  switch (inAxis)
60  {
61  case Axis::ROW:
62  {
63  const auto inShape = inArray.shape();
64  if (inShape.rows < 2)
65  {
66  THROW_INVALID_ARGUMENT_ERROR("input array must have more than 1 row.");
67  }
68 
69  // first do the first and last rows
70  auto returnArray = NdArray<double>(inShape);
71  for (uint32 col = 0; col < inShape.cols; ++col)
72  {
73  returnArray(0, col) = static_cast<double>(inArray(1, col)) - static_cast<double>(inArray(0, col));
74  returnArray(-1, col) = static_cast<double>(inArray(-1, col)) - static_cast<double>(inArray(-2, col));
75  }
76 
77  // then rip through the rest of the array
78  for (uint32 col = 0; col < inShape.cols; ++col)
79  {
80  for (uint32 row = 1; row < inShape.rows - 1; ++row)
81  {
82  returnArray(row, col) = (static_cast<double>(inArray(row + 1, col)) - static_cast<double>(inArray(row - 1, col))) / 2.0;
83  }
84  }
85 
86  return returnArray;
87  }
88  case Axis::COL:
89  {
90  const auto inShape = inArray.shape();
91  if (inShape.cols < 2)
92  {
93  THROW_INVALID_ARGUMENT_ERROR("input array must have more than 1 columns.");
94  }
95 
96  // first do the first and last columns
97  auto returnArray = NdArray<double>(inShape);
98  for (uint32 row = 0; row < inShape.rows; ++row)
99  {
100  returnArray(row, 0) = static_cast<double>(inArray(row, 1)) - static_cast<double>(inArray(row, 0));
101  returnArray(row, -1) = static_cast<double>(inArray(row, -1)) - static_cast<double>(inArray(row, -2));
102  }
103 
104  // then rip through the rest of the array
105  for (uint32 row = 0; row < inShape.rows; ++row)
106  {
107  for (uint32 col = 1; col < inShape.cols - 1; ++col)
108  {
109  returnArray(row, col) = (static_cast<double>(inArray(row, col + 1)) - static_cast<double>(inArray(row, col - 1))) / 2.0;
110  }
111  }
112 
113  return returnArray;
114  }
115  default:
116  {
117  // will return the gradient of the flattened array
118  if (inArray.size() < 2)
119  {
120  THROW_INVALID_ARGUMENT_ERROR("input array must have more than 1 element.");
121  }
122 
123  auto returnArray = NdArray<double>(1, inArray.size());
124  returnArray[0] = static_cast<double>(inArray[1]) - static_cast<double>(inArray[0]);
125  returnArray[-1] = static_cast<double>(inArray[-1]) - static_cast<double>(inArray[-2]);
126 
127  stl_algorithms::transform(inArray.cbegin() + 2, inArray.cend(), inArray.cbegin(), returnArray.begin() + 1,
128  [](dtype value1, dtype value2) -> double
129  {
130  return (static_cast<double>(value1) - static_cast<double>(value2)) / 2.0;
131  });
132 
133  return returnArray;
134  }
135  }
136  }
137 
138  //============================================================================
139  // Method Description:
149  template<typename dtype>
150  NdArray<std::complex<double>> gradient(const NdArray<std::complex<dtype>>& inArray, Axis inAxis = Axis::ROW)
151  {
153 
154  switch (inAxis)
155  {
156  case Axis::ROW:
157  {
158  const auto inShape = inArray.shape();
159  if (inShape.rows < 2)
160  {
161  THROW_INVALID_ARGUMENT_ERROR("input array must have more than 1 row.");
162  }
163 
164  // first do the first and last rows
165  auto returnArray = NdArray<std::complex<double>>(inShape);
166  for (uint32 col = 0; col < inShape.cols; ++col)
167  {
168  returnArray(0, col) = complex_cast<double>(inArray(1, col)) - complex_cast<double>(inArray(0, col));
169  returnArray(-1, col) = complex_cast<double>(inArray(-1, col)) - complex_cast<double>(inArray(-2, col));
170  }
171 
172  // then rip through the rest of the array
173  for (uint32 col = 0; col < inShape.cols; ++col)
174  {
175  for (uint32 row = 1; row < inShape.rows - 1; ++row)
176  {
177  returnArray(row, col) = (complex_cast<double>(inArray(row + 1, col)) -
178  complex_cast<double>(inArray(row - 1, col))) / 2.0;
179  }
180  }
181 
182  return returnArray;
183  }
184  case Axis::COL:
185  {
186  const auto inShape = inArray.shape();
187  if (inShape.cols < 2)
188  {
189  THROW_INVALID_ARGUMENT_ERROR("input array must have more than 1 columns.");
190  }
191 
192  // first do the first and last columns
193  auto returnArray = NdArray<std::complex<double>>(inShape);
194  for (uint32 row = 0; row < inShape.rows; ++row)
195  {
196  returnArray(row, 0) = complex_cast<double>(inArray(row, 1)) - complex_cast<double>(inArray(row, 0));
197  returnArray(row, -1) = complex_cast<double>(inArray(row, -1)) - complex_cast<double>(inArray(row, -2));
198  }
199 
200  // then rip through the rest of the array
201  for (uint32 row = 0; row < inShape.rows; ++row)
202  {
203  for (uint32 col = 1; col < inShape.cols - 1; ++col)
204  {
205  returnArray(row, col) = (complex_cast<double>(inArray(row, col + 1)) -
206  complex_cast<double>(inArray(row, col - 1))) / 2.0;
207  }
208  }
209 
210  return returnArray;
211  }
212  default:
213  {
214  // will return the gradient of the flattened array
215  if (inArray.size() < 2)
216  {
217  THROW_INVALID_ARGUMENT_ERROR("input array must have more than 1 element.");
218  }
219 
220  auto returnArray = NdArray<std::complex<double>>(1, inArray.size());
221  returnArray[0] = complex_cast<double>(inArray[1]) - complex_cast<double>(inArray[0]);
222  returnArray[-1] = complex_cast<double>(inArray[-1]) - complex_cast<double>(inArray[-2]);
223 
224  stl_algorithms::transform(inArray.cbegin() + 2, inArray.cend(), inArray.cbegin(), returnArray.begin() + 1,
225  [](const std::complex<dtype>& value1, const std::complex<dtype>& value2) -> std::complex<double>
226  {
227  return (complex_cast<double>(value1) - complex_cast<double>(value2)) / 2.0;
228  });
229 
230  return returnArray;
231  }
232  }
233  }
234 } // namespace nc
#define THROW_INVALID_ARGUMENT_ERROR(msg)
Definition: Error.hpp:36
#define STATIC_ASSERT_ARITHMETIC(dtype)
Definition: StaticAsserts.hpp:37
size_type size() const noexcept
Definition: NdArrayCore.hpp:4296
const_iterator cbegin() const noexcept
Definition: NdArrayCore.hpp:1216
Shape shape() const noexcept
Definition: NdArrayCore.hpp:4283
const_iterator cend() const noexcept
Definition: NdArrayCore.hpp:1524
iterator begin() noexcept
Definition: NdArrayCore.hpp:1166
OutputIt transform(InputIt first, InputIt last, OutputIt destination, UnaryOperation unaryFunction)
Definition: StlAlgorithms.hpp:702
Definition: Coordinate.hpp:45
Axis
Enum To describe an axis.
Definition: Types.hpp:46
NdArray< double > gradient(const NdArray< dtype > &inArray, Axis inAxis=Axis::ROW)
Definition: gradient.hpp:55
std::uint32_t uint32
Definition: Types.hpp:40