MaskedArray

class numpy::MaskedArray

numpy C++ class.

Example

#include <numpy/np_ndarray.h>
using namespace numpy;

// Use MaskedArray
MaskedArray obj;
// ... operations ...

Constructors

Signature	Location	Example
MaskedArray(const NDArray<T>&data, Tfill_value = get_default_fill_value<T>())	NP_MASKED_ARRAY.H:54
MaskedArray(const NDArray<T>&data, const NDArray<bool\_>&mask, Tfill_value = get_default_fill_value<T>())	NP_MASKED_ARRAY.H:60
MaskedArray(const MaskedArray<T>&other)	NP_MASKED_ARRAY.H:69

Operators

Signature	Return Type	Location	Example
MaskedArray<T>& operator=(const MaskedArray<T>&other)	MaskedArray<T>&	NP_MASKED_ARRAY.H:73
MaskedArray<T> operator+(const MaskedArray<T>&other)	MaskedArray<T>	NP_MASKED_ARRAY.H:315
MaskedArray<T> operator+(Tscalar)	MaskedArray<T>	NP_MASKED_ARRAY.H:348
MaskedArray<T> operator-(const MaskedArray<T>&other)	MaskedArray<T>	NP_MASKED_ARRAY.H:362
MaskedArray<T> operator-(Tscalar)	MaskedArray<T>	NP_MASKED_ARRAY.H:394
MaskedArray<T> operator\*(const MaskedArray<T>&other)	MaskedArray<T>	NP_MASKED_ARRAY.H:408
MaskedArray<T> operator\*(Tscalar)	MaskedArray<T>	NP_MASKED_ARRAY.H:440
MaskedArray<T> operator/(const MaskedArray<T>&other)	MaskedArray<T>	NP_MASKED_ARRAY.H:454
MaskedArray<T> operator/(Tscalar)	MaskedArray<T>	NP_MASKED_ARRAY.H:500
MaskedArray<T> operator-()	MaskedArray<T>	NP_MASKED_ARRAY.H:538
MaskedArray<bool\_> operator==(const MaskedArray<T>&other)	MaskedArray<bool_>	NP_MASKED_ARRAY.H:587
MaskedArray<bool\_> operator<(const MaskedArray<T>&other)	MaskedArray<bool_>	NP_MASKED_ARRAY.H:618
MaskedArray<bool\_> operator>(const MaskedArray<T>&other)	MaskedArray<bool_>	NP_MASKED_ARRAY.H:649
MaskedArray<bool\_> operator!=(const MaskedArray<T>&other)	MaskedArray<bool_>	NP_MASKED_ARRAY.H:680
MaskedArray<bool\_> operator<=(const MaskedArray<T>&other)	MaskedArray<bool_>	NP_MASKED_ARRAY.H:711
MaskedArray<bool\_> operator>=(const MaskedArray<T>&other)	MaskedArray<bool_>	NP_MASKED_ARRAY.H:742

Indexing / Selection

Signature	Return Type	Location	Example
T getElementAt(const std::vector<size_t>&indices)	T	NP_MASKED_ARRAY.H:113	View
size_t getSize()	size_t	NP_MASKED_ARRAY.H:106	View
U get_default_fill_value()	U	NP_MASKED_ARRAY.H:151

Statistics

Signature	Return Type	Location	Example
T max()	T	NP_MASKED_ARRAY.H:817	View
double mean()	double	NP_MASKED_ARRAY.H:791	View
T min()	T	NP_MASKED_ARRAY.H:800	View
double std()	double	NP_MASKED_ARRAY.H:834	View
T sum()	T	NP_MASKED_ARRAY.H:777	View
double var()	double	NP_MASKED_ARRAY.H:857	View

Type Checking

Signature	Return Type	Location	Example
bool is_all_masked()	bool	NP_MASKED_ARRAY.H:174
bool is_masked()	bool	NP_MASKED_ARRAY.H:163	View
bool is_masked_at(const std::vector<size_t>&indices)	bool	NP_MASKED_ARRAY.H:131	View

Other Methods

Signature	Return Type	Location	Example
size_t count()	size_t	NP_MASKED_ARRAY.H:185	View
const NDArray<T>& data()	const NDArray<T>&	NP_MASKED_ARRAY.H:88	View
NDArray<T>& data()	NDArray<T>&	NP_MASKED_ARRAY.H:89	View
T fill_value()	T	NP_MASKED_ARRAY.H:96	View
bool hard_mask()	bool	NP_MASKED_ARRAY.H:100
void harden_mask()	void	NP_MASKED_ARRAY.H:101	View
const NDArray<bool\_>& mask()	const NDArray<bool_>&	NP_MASKED_ARRAY.H:92	View
NDArray<bool\_>& mask()	NDArray<bool_>&	NP_MASKED_ARRAY.H:93	View
void mask_at(const std::vector<size_t>&indices)	void	NP_MASKED_ARRAY.H:136	View
void mask_equal(Tvalue)	void	NP_MASKED_ARRAY.H:289
void mask_greater(Tvalue)	void	NP_MASKED_ARRAY.H:294
void mask_invalid()	void	NP_MASKED_ARRAY.H:280
void mask_less(Tvalue)	void	NP_MASKED_ARRAY.H:299
void mask_where(Funccondition)	void	NP_MASKED_ARRAY.H:247	View
void mask_where(const NDArray<bool\_>&condition)	void	NP_MASKED_ARRAY.H:262	View
T prod()	T	NP_MASKED_ARRAY.H:863	View
void setElementAt(const std::vector<size_t>&indices, Tvalue)	void	NP_MASKED_ARRAY.H:121	View
void set_fill_value(Tvalue)	void	NP_MASKED_ARRAY.H:97	View
void soften_mask()	void	NP_MASKED_ARRAY.H:102	View
void unmask_all()	void	NP_MASKED_ARRAY.H:304
void unmask_at(const std::vector<size_t>&indices)	void	NP_MASKED_ARRAY.H:140	View

Code Examples

The following examples are extracted from the test suite.

getElementAt (np_test_1_all.cpp:49)

using namespace numpy::benchmark;
using namespace numpy::char_;

// Helper functions for array comparison tests
namespace {
    const double TOLERANCE = 1e-10;

    bool allTrue(const NDArray<bool>& arr) {
        std::vector<size_t> indices(arr.getShape().size(), 0);
        do {
            if (!arr.getElementAt(indices)) {
                return false;
            }
        } while (incrementIndices(indices, arr.getShape()));
        return true;
    }

    bool allFalse(const NDArray<bool>& arr) {
        std::vector<size_t> indices(arr.getShape().size(), 0);
        do {
            if (arr.getElementAt(indices)) {

getSize (np_test_1_all.cpp:2172)

    }

    std::cout << " -> tests passed" << std::endl;
}

void testArrayEdgeCases() {
    std::cout << "========= testArrayEdgeCases =======================";

    // Test 1: Empty arrays
    NDArray<object_> empty_array = createObjectZeros({0});
    if (!(empty_array.getSize() == 0)) {
        std::string description = std::string("testArrayEdgeCases():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(empty_array.getSize() == 0)";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    if (!(empty_array.getShape() == std::vector<size_t>{0})) {
        std::string description = std::string("testArrayEdgeCases():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(empty_array.getShape() == std::vector<size_t>{0})";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }

max (np_test_1_all.cpp:7274)

    if (sizeof(uintp) == sizeof(void*)) {
        // std::cout << "                -> uintp size matches pointer size";
    } else {
        // std::cout << "  ✗ uintp size doesn't match pointer size" << std::endl;
    }

    // Test range limits
    // std::cout << "Range Information:" << std::endl;
    // std::cout << "  intp min: " << std::numeric_limits<intp>::min() << std::endl;
    // std::cout << "  intp max: " << std::numeric_limits<intp>::max() << std::endl;
    // std::cout << "  uintp max: " << std::numeric_limits<uintp>::max() << std::endl;
    // std::cout << "  longdouble digits: " << std::numeric_limits<longdouble>::digits << std::endl;

    std::cout << " -> tests passed" << std::endl;
}

void testComplexArithmeticExtendedTypes() {
    std::cout << "========= testComplexArithmeticExtendedTypes =======================";

    clongdouble c1(3.0L, 4.0L);  // 3 + 4i

mean (np_test_1_all.cpp:11714)

    // Create test array
    auto array = createInt32Array({ 2, 3 }, 0);
    array.setElementAt({ 0, 0 }, 1);
    array.setElementAt({ 0, 1 }, 2);
    array.setElementAt({ 0, 2 }, 3);
    array.setElementAt({ 1, 0 }, 4);
    array.setElementAt({ 1, 1 }, 5);
    array.setElementAt({ 1, 2 }, 6);

    // Test mean without axis
    auto mean_all = mean(array);
    if (!(approx_equal(mean_all.getElementAt({ 0 }), 3.5, 1e-10))) {
        std::string description = std::string("testBasicStatistics():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(approx_equal(mean_all.getElementAt({ 0 }), 3.5, 1e-10))";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    // std::cout << "[OK] mean (all elements) works correctly\n";

    // Test mean along axis 0
    auto mean_axis0 = mean(array, 0);
    if (!(mean_axis0.getShape()[0] == 3)) {

min (np_test_1_all.cpp:2350)

        if (i % 3 == 0) {
            large_array.setElementAt({i}, object_(static_cast<int>(i)));
        } else if (i % 3 == 1) {
            large_array.setElementAt({i}, object_(static_cast<double>(i) * 0.5));
        } else {
            large_array.setElementAt({i}, object_(std::string("str") + std::to_string(i)));
        }
    }

    // Verify pattern
    for (size_t i = 0; i < std::min(large_size, size_t(100)); ++i) {  // Check first 100
        object_ obj = large_array.getElementAt({i});
        if (i % 3 == 0) {
            if (!(obj.is_type<int>())) {
                std::string description = std::string("testArrayEdgeCases():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(obj.is_type<int>())";
                std::cout << std::string("[FAIL] ") + description << std::endl;
                throw std::runtime_error(description);
            }
        } else if (i % 3 == 1) {
            if (!(obj.is_type<double>())) {
                std::string description = std::string("unknown_function():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(obj.is_type<double>())";

std (np_test_1_all.cpp:11836)

    auto var_sample = var(array, std::nullopt, 1);  // ddof=1 for sample variance
    double expected_sample_var = 1.25 * 4.0 / 3.0;  // Bessel's correction
    if (!(approx_equal(var_sample.getElementAt({ 0 }), expected_sample_var, 1e-10))) {
        std::string description = std::string("testVarianceAndStandardDeviation():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(approx_equal(var_sample.getElementAt({ 0 }), expected_sample_var, 1e-10))";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    // std::cout << "[OK] Sample variance works correctly\n";

    // Test standard deviation
    auto std_result = numpy::std(array, std::nullopt, 0);
    if (!(approx_equal(std_result.getElementAt({ 0 }), std::sqrt(expected_var), 1e-10))) {
        std::string description = std::string("testVarianceAndStandardDeviation():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(approx_equal(std_result.getElementAt({ 0 }), std::sqrt(expected_var), 1e-10))";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    // std::cout << "[OK] Standard deviation works correctly\n";

    // Test along axis
    auto array2d = createFloat64Array({ 2, 2 }, 0);
    array2d.setElementAt({ 0, 0 }, 1.0);

sum (np_test_1_all.cpp:11766)

        throw std::runtime_error(description);
    }
    if (!(approx_equal(mean_axis1.getElementAt({ 1 }), 5.0, 1e-10))) {
        std::string description = std::string("testBasicStatistics():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(approx_equal(mean_axis1.getElementAt({ 1 }), 5.0, 1e-10))";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    // std::cout << "[OK] mean along axis 1 works correctly\n";

    // Test sum
    auto sum_all = sum(array);
    if (!(sum_all.getElementAt({ 0 }) == 21)) {
        std::string description = std::string("testBasicStatistics():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(sum_all.getElementAt({ 0 }) == 21)";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    // std::cout << "[OK] sum works correctly\n";

    // Test min and max
    auto min_all = min(array);
    auto max_all = max(array);

var (np_test_1_all.cpp:11816)

    std::cout << "========= testVarianceAndStandardDeviation =======================";

    // Create test array with known variance
    auto array = createFloat64Array({ 4 }, 0);
    array.setElementAt({ 0 }, 1.0);
    array.setElementAt({ 1 }, 2.0);
    array.setElementAt({ 2 }, 3.0);
    array.setElementAt({ 3 }, 4.0);

    // Test variance (population)
    auto var_result = var(array, std::nullopt, 0);  // ddof=0 for population variance
    double expected_var = 1.25;  // Known variance for [1,2,3,4]
    if (!(approx_equal(var_result.getElementAt({ 0 }), expected_var, 1e-10))) {
        std::string description = std::string("testVarianceAndStandardDeviation():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(approx_equal(var_result.getElementAt({ 0 }), expected_var, 1e-10))";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    // std::cout << "[OK] Population variance works correctly\n";

    // Test variance (sample)
    auto var_sample = var(array, std::nullopt, 1);  // ddof=1 for sample variance

is_masked (np_test_5_all.cpp:3746)

      std::cout << "========= is_masked =======================";

      // Test with masked array that has some masked values
      numpy::NDArray<double> arr1({ 3 });
      arr1.setElementAt({ 0 }, 1.0);
      arr1.setElementAt({ 1 }, std::numeric_limits<double>::quiet_NaN());
      arr1.setElementAt({ 2 }, 3.0);

      auto ma1 = numpy::masked_invalid(arr1);

      if (!numpy::is_masked(ma1)) {
        std::cout << "  [FAIL] : in np_test_is_masked() : Array with NaN not detected as masked";
        throw std::runtime_error("Test failed");
      }

      // Test with masked array with no masked values
      numpy::NDArray<double> arr2({ 3 });
      arr2.setElementAt({ 0 }, 1.0);
      arr2.setElementAt({ 1 }, 2.0);
      arr2.setElementAt({ 2 }, 3.0);

is_masked_at (np_test_2_all.cpp:6004)

      data.setElementAt({ 2, 0 }, 7.0);
      data.setElementAt({ 2, 1 }, 8.0);
      data.setElementAt({ 2, 2 }, 9.0);

      auto mask = createBoolArray({ 3, 3 }, false);
      mask.setElementAt({ 0, 0 }, true);
      mask.setElementAt({ 1, 1 }, true);

      MaskedArray<double> ma(data, mask, -999.0);

      if (!(ma.is_masked_at({ 0, 0 }))) {
          std::string description = std::string("testMaskedArrayCreation():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(ma.is_masked_at({ 0, 0 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }
      if (!(ma.is_masked_at({ 1, 1 }))) {
          std::string description = std::string("testMaskedArrayCreation():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(ma.is_masked_at({ 1, 1 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }
      if (!(!ma.is_masked_at({ 0, 1 }))) {

count (np_test_1_all.cpp:3616)

    // Create larger arrays for performance testing
    auto large_arr = NDArray<double>::createOnes({100, 100});
    auto broadcast_arr = NDArray<double>::createOnes({1, 100});

    // Time the operation (basic timing)
    auto start = std::chrono::high_resolution_clock::now();
    auto result = large_arr.addArrays(broadcast_arr);
    auto end = std::chrono::high_resolution_clock::now();

    auto duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start);
    // std::cout << "Large array broadcasting took " << duration.count() << " microseconds" << std::endl;

    // Verify result shape
    if (!((result.getShape() == std::vector<size_t>{100, 100}))) {
        std::string description = std::string("test_broadcasting_performance():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !((result.getShape() == std::vector<size_t>{100, 100}))";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    if (!(result.getElementAt({50, 50}) == 2.0)) {
        std::string description = std::string("test_broadcasting_performance():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(result.getElementAt({50, 50}) == 2.0)";
        std::cout << std::string("[FAIL] ") + description << std::endl;

data (np_test_1_all.cpp:2084)

    }
    if (!(derived_ptr->extra == 2)) {
        std::string description = std::string("testObjectHolderEdgeCases():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(derived_ptr->extra == 2)";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }

    // Test 6: Large object handling
    struct LargeObject {
        std::vector<double> data;
        LargeObject() : data(10000, 3.14159) {}  // Large object
    };

    LargeObject large;
    object_ large_obj(large);  // Should handle large objects
    if (!(!large_obj.is_null())) {
        std::string description = std::string("testObjectHolderEdgeCases():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(!large_obj.is_null())";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    if (!(large_obj.is_type<LargeObject>())) {

data (np_test_1_all.cpp:2084)

    }
    if (!(derived_ptr->extra == 2)) {
        std::string description = std::string("testObjectHolderEdgeCases():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(derived_ptr->extra == 2)";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }

    // Test 6: Large object handling
    struct LargeObject {
        std::vector<double> data;
        LargeObject() : data(10000, 3.14159) {}  // Large object
    };

    LargeObject large;
    object_ large_obj(large);  // Should handle large objects
    if (!(!large_obj.is_null())) {
        std::string description = std::string("testObjectHolderEdgeCases():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(!large_obj.is_null())";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    if (!(large_obj.is_type<LargeObject>())) {

fill_value (np_test_3_all.cpp:21316)

    std::cout << "========= Pad jagged array ====" ;

    // Create jagged array
    numpy::NDArray<numpy::object_> jagged = numpy::make_jagged_array({
        {1, 2, 3},
        {4, 5},
        {6}
    });

    // Pad with -1
    numpy::object_ fill_value(-1);
    numpy::NDArray<numpy::object_> padded = numpy::NestedArray<int>::padJagged(jagged, fill_value);

    // Should be 3x3
    CHECK(padded.getShape().size() == 2, "np_test_nested_pad", "ndim != 2");
    CHECK(padded.getShape()[0] == 3, "np_test_nested_pad", "shape[0] != 3");
    CHECK(padded.getShape()[1] == 3, "np_test_nested_pad", "shape[1] != 3");

    // Check values - must store object_ to keep pointers valid (avoid dangling)
    numpy::object_ obj00 = padded.getElementAt({0, 0});
    numpy::object_ obj01 = padded.getElementAt({0, 1});

harden_mask (np_test_2_all.cpp:6446)

          std::string description = std::string("testHardSoftMasks():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(!ma_soft.is_masked_at({ 1 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }
      // std::cout << "[OK] Soft mask allows unmasking\n";

      // Test hard mask - need to set hard mask mode first
      MaskedArray<double> ma_hard(data, mask, -999.0);

      // std::cout << ma_hard.is_masked_at({ 1 }) << std::endl;
      ma_hard.harden_mask();
      ma_hard.unmask_at({ 1 });
      if (!(ma_hard.is_masked_at({ 1 }))) {
          std::string description = std::string("testHardSoftMasks():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(ma_hard.is_masked_at({ 1 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }
      // std::cout << "[OK] Hard mask prevents unmasking\n";

      std::cout << " -> tests passed\n";
    }

mask (np_test_1_all.cpp:27691)

    void test_mask_rows_cols() {
      std::cout << "========= ma::mask_rows() and ma::mask_cols(): mask operations ===";

      auto data = numpy::array({ {1.0, 2.0, 3.0}, {4.0, 5.0, 6.0} });
      auto mask = numpy::array({ {false, true, false}, {false, false, false} });
      auto ma = numpy::ma::masked_array(data, mask);

      auto result_rows = numpy::ma::mask_rows(ma);
      auto result_cols = numpy::ma::mask_cols(ma);

      // std::cout << "  ma::mask_rows(): row 0 fully masked: " << (result_rows.mask().getElementAt({ 0, 0 }) ? "true" : "false") << std::endl;
      // std::cout << "  ma::mask_cols(): column 1 fully masked: " << (result_cols.mask().getElementAt({ 1, 1 }) ? "true" : "false");
      std::cout << " -> tests passed";
    }

    void test_compress_rowcols() {
      std::cout << "========= ma::compress_rowcols(): remove masked rows/columns ===";

      auto data = numpy::array({ {1.0, 2.0, 3.0}, {4.0, 5.0, 6.0}, {7.0, 8.0, 9.0} });
      auto mask = numpy::array({ {false, false, false}, {false, true, false}, {false, false, false} });
      auto ma = numpy::ma::masked_array(data, mask);

mask (np_test_1_all.cpp:27691)

    void test_mask_rows_cols() {
      std::cout << "========= ma::mask_rows() and ma::mask_cols(): mask operations ===";

      auto data = numpy::array({ {1.0, 2.0, 3.0}, {4.0, 5.0, 6.0} });
      auto mask = numpy::array({ {false, true, false}, {false, false, false} });
      auto ma = numpy::ma::masked_array(data, mask);

      auto result_rows = numpy::ma::mask_rows(ma);
      auto result_cols = numpy::ma::mask_cols(ma);

      // std::cout << "  ma::mask_rows(): row 0 fully masked: " << (result_rows.mask().getElementAt({ 0, 0 }) ? "true" : "false") << std::endl;
      // std::cout << "  ma::mask_cols(): column 1 fully masked: " << (result_cols.mask().getElementAt({ 1, 1 }) ? "true" : "false");
      std::cout << " -> tests passed";
    }

    void test_compress_rowcols() {
      std::cout << "========= ma::compress_rowcols(): remove masked rows/columns ===";

      auto data = numpy::array({ {1.0, 2.0, 3.0}, {4.0, 5.0, 6.0}, {7.0, 8.0, 9.0} });
      auto mask = numpy::array({ {false, false, false}, {false, true, false}, {false, false, false} });
      auto ma = numpy::ma::masked_array(data, mask);

mask_at (np_test_5_all.cpp:4074)

      // Create a 2D array
      numpy::NDArray<double> arr({ 3, 4 });
      for (size_t i = 0; i < 3; ++i) {
        for (size_t j = 0; j < 4; ++j) {
          arr.setElementAt({ i, j }, static_cast<double>(i * 4 + j));
        }
      }
      numpy::MaskedArray<double> ma(arr);

      // Mask one element at (0, 1)
      ma.mask_at({ 0, 1 });

      // Mask entire row containing masked element (axis=0)
      auto ma_rows = numpy::mask_rowcols(ma, 0);

      // Check that entire first row is masked
      for (size_t j = 0; j < 4; ++j) {
        if (!ma_rows.mask().getElementAt({ 0, j })) {
          std::cout << "  [FAIL] : in np_test_mask_rowcols() : Row 0 should be fully masked";
          throw std::runtime_error("Test failed");
        }

mask_where (np_test_2_all.cpp:6308)

      data.setElementAt({ 3 }, 40.0);

      auto mask = createBoolArray({ 4 }, false);
      MaskedArray<double> ma(data, mask, -999.0);

      // Test mask_where
      auto condition = createBoolArray({ 4 }, false);
      condition.setElementAt({ 0 }, true);
      condition.setElementAt({ 2 }, true);

      ma.mask_where(condition);
      if (!(ma.is_masked_at({ 0 }))) {
          std::string description = std::string("testMaskOperations():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(ma.is_masked_at({ 0 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }
      if (!(!ma.is_masked_at({ 1 }))) {
          std::string description = std::string("testMaskOperations():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(!ma.is_masked_at({ 1 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }

mask_where (np_test_2_all.cpp:6308)

      data.setElementAt({ 3 }, 40.0);

      auto mask = createBoolArray({ 4 }, false);
      MaskedArray<double> ma(data, mask, -999.0);

      // Test mask_where
      auto condition = createBoolArray({ 4 }, false);
      condition.setElementAt({ 0 }, true);
      condition.setElementAt({ 2 }, true);

      ma.mask_where(condition);
      if (!(ma.is_masked_at({ 0 }))) {
          std::string description = std::string("testMaskOperations():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(ma.is_masked_at({ 0 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }
      if (!(!ma.is_masked_at({ 1 }))) {
          std::string description = std::string("testMaskOperations():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(!ma.is_masked_at({ 1 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }

prod (np_test_1_all.cpp:11795)

        throw std::runtime_error(description);
    }
    // std::cout << "[OK] min and max work correctly\n";

    // Test product
    auto small_array = createInt32Array({ 3 }, 0);
    small_array.setElementAt({ 0 }, 2);
    small_array.setElementAt({ 1 }, 3);
    small_array.setElementAt({ 2 }, 4);

    auto prod_result = prod(small_array);
    if (!(prod_result.getElementAt({ 0 }) == 24)) {
        std::string description = std::string("testBasicStatistics():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(prod_result.getElementAt({ 0 }) == 24)";
        std::cout << std::string("[FAIL] ") + description << std::endl;
        throw std::runtime_error(description);
    }
    // std::cout << "[OK] prod works correctly\n";
    std::cout << " -> tests passed" << std::endl;
  }

  void testVarianceAndStandardDeviation() {

setElementAt (np_test_1_all.cpp:135)

    //ones_array.printArray();

    std::cout << " -> tests passed" << std::endl;
}

void testElementAccess() {
    std::cout << "========= testElementAccess =======================";

    auto array = createInt32Array({3, 3}, 0);

    array.setElementAt({0, 0}, 1);
    array.setElementAt({0, 1}, 2);
    array.setElementAt({0, 2}, 3);
    array.setElementAt({1, 1}, 5);
    array.setElementAt({2, 2}, 9);

    // std::cout << "Array after setting elements:" << std::endl;
    //array.printArray();

    // std::cout << "Element at (1,1): " << array.getElementAt({1, 1}) << std::endl;
    // std::cout << "Element at (2,2): " << array.getElementAt({2, 2});

set_fill_value (np_test_5_all.cpp:3829)

    void np_test_common_fill_value() {
      std::cout << "========= common_fill_value =======================";

      // Create two masked arrays with the same fill value
      numpy::NDArray<double> arr1({ 3 });
      arr1.setElementAt({ 0 }, 1.0);
      arr1.setElementAt({ 1 }, 2.0);
      arr1.setElementAt({ 2 }, 3.0);
      numpy::MaskedArray<double> ma1(arr1);
      ma1.set_fill_value(999.0);

      numpy::NDArray<double> arr2({ 3 });
      arr2.setElementAt({ 0 }, 4.0);
      arr2.setElementAt({ 1 }, 5.0);
      arr2.setElementAt({ 2 }, 6.0);
      numpy::MaskedArray<double> ma2(arr2);
      ma2.set_fill_value(999.0);

      auto fill = numpy::common_fill_value(ma1, ma2);

soften_mask (np_test_5_all.cpp:4116)

      for (size_t i = 0; i < 5; ++i) {
        arr.setElementAt({ i }, static_cast<double>(i + 1));
      }
      numpy::MaskedArray<double> ma(arr);

      // Mask some elements
      ma.mask_at({ 1 });
      ma.mask_at({ 3 });

      // Soften the mask (allows modification)
      ma.soften_mask();

      // After softening, we should still have the masked elements
      if (!ma.mask().getElementAt({ 1 }) || !ma.mask().getElementAt({ 3 })) {
        std::cout << "  [FAIL] : in np_test_soften_mask() : Masked elements lost after soften";
        throw std::runtime_error("Test failed");
      }

      // Verify unmasked elements are still unmasked
      if (ma.mask().getElementAt({ 0 }) || ma.mask().getElementAt({ 2 }) || ma.mask().getElementAt({ 4 })) {
        std::cout << "  [FAIL] : in np_test_soften_mask() : Unmasked elements incorrectly masked";

unmask_at (np_test_2_all.cpp:6335)

          std::string description = std::string("testMaskOperations():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(!ma.is_masked_at({ 3 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }
      // std::cout << "[OK] mask_where() applies conditional masking\n";

      // Test unmask (for soft masks)
      auto unmask_condition = createBoolArray({ 4 }, false);
      unmask_condition.setElementAt({ 0 }, true);

      ma.unmask_at({ 0 });
      if (!(!ma.is_masked_at({ 0 }))) {
          std::string description = std::string("testMaskOperations():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(!ma.is_masked_at({ 0 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }
      if (!(ma.is_masked_at({ 2 }))) {
          std::string description = std::string("testMaskOperations():") + __FILE__ + ":" + std::to_string(__LINE__) + ": !(ma.is_masked_at({ 2 }))";
          std::cout << std::string("[FAIL] ") + description << std::endl;
          throw std::runtime_error(description);
      }