初学C编写Python C扩展模块返回inf/极大值问题求助

我刚开始学习C语言，这是我用C编码的第一天，还请多多包涵。我尝试编写一个Python模块，接收numpy数组，对其进行滑动窗口处理，为每个子窗口计算标准差，并将结果存储到新的numpy数组中返回给Python。不知为何，该模块在Python中返回的是'inf'或极大数值（具体结果随代码的随机修改而变化）。我尝试将计算简化为仅求均值（相当于一个简易低通滤波器），但仍得到这些不合理的数值。我怀疑是在与Python进行值传递和返回的环节中存在疏漏。

C代码 #

#include <Python.h>
#include <stdio.h>
#include <stdlib.h>
#include <unistd.h>
#include <numpy/arrayobject.h>


float calculateStandardDeviation(float *arr, int startIndex, int windowSize) {
    float mean = 0.0;
    float variance = 0.0;

    // Calculate the mean of the subarray
    for (int i = startIndex; i < startIndex + windowSize; i++) {
        mean += arr[i];
    }
    mean /= windowSize;

    // Calculate the variance
    for (int i = startIndex; i < startIndex + windowSize; i++) {
        variance += pow(arr[i] - mean, 2);
    }
    variance /= windowSize;

    // Calculate the standard deviation
    float standardDeviation = sqrt(variance);
    return standardDeviation;
}

static PyArrayObject *method_stdarray(PyObject * self, PyObject *args){
    PyArrayObject *inputArrayObj;
    int windowSize, stepSize;

    if (!PyArg_ParseTuple(args, "Oii", &inputArrayObj, &windowSize, &stepSize)) {
        return NULL;
    }
    //Make sure it's a numpy array
    if (!PyArray_Check(inputArrayObj)) {
        PyErr_SetString(PyExc_TypeError, "Input array must be numpy array.");
        return NULL;
    }
    // Get dims of array
    int ndim = PyArray_NDIM(inputArrayObj);

    // Make sure its a 1D-array
    if (ndim > 1){PyErr_SetString(PyExc_BaseException, "Input numpy array has too many dimensions."); return NULL;}

    if (PyArray_TYPE(inputArrayObj) != NPY_DOUBLE) {
        PyErr_SetString(PyExc_ValueError, "Input array must be of type float64.");
        return NULL;
    }

    //Get the size
    int arraySize = PyArray_Size(inputArrayObj);

    npy_intp size = PyArray_SIZE(inputArrayObj);
    double* inputArray = PyArray_DATA(inputArrayObj);


    if (inputArray == NULL) {
        PyErr_NoMemory();
        return NULL;
    }


    int numWindows = (arraySize - windowSize) / stepSize + 1;
    npy_intp numpyWindows = numWindows;
    //PyObject *resultList = PyList_New(numWindows);

    PyArrayObject* resultArray = (PyArrayObject*)PyArray_EMPTY(1, &numpyWindows, NPY_DOUBLE, 0);

    if (resultArray == NULL) {
        PyErr_SetString(PyExc_MemoryError, "Failed to allocate memory for the result array.");
        return NULL;
    }

    double* resultData = PyArray_DATA(resultArray);
    
    for (int i = 0; i < numWindows; i++) {
        int startIndex = i * stepSize;
        float result = calculateStandardDeviation(inputArray, startIndex, windowSize);
        resultData[i] = result;
    }
    
    Py_INCREF(inputArrayObj);

    return (PyObject*)resultArray;

}


static PyMethodDef meanyMethods[] = {
    {"vararray", method_stdarray, METH_VARARGS, "Takes a list and two integers, returns a new list of variances of subwindows"},
    {NULL, NULL, 0, NULL}
};


static struct PyModuleDef meanymodule = {
    PyModuleDef_HEAD_INIT,
    "meany",
    "My set of functions for analysing mea data",
    -1,
    meanyMethods
};

PyMODINIT_FUNC PyInit_meany(void) {
    import_array();
    return PyModule_Create(&meanymodule);
}

Python调用代码 #

import meany
import random
import numpy as np

fakeList = np.ones(100)
for i in range(0,100):
    fakeList[i] = fakeList[i] + (random.random())
stdList = meany.vararray(fakeList,20,1)

print(stdList)

问题分析与修复方案 #

核心问题 #

1. 类型不匹配：numpy数组指定为NPY_DOUBLE类型（对应C的double），但calculateStandardDeviation函数接收的是float *arr参数。将double*强制转换为float*会导致内存解析错误，读取到完全错误的数值，最终计算出inf或极大值。
2. 缺少头文件：代码中使用了pow和sqrt函数，但未包含<math.h>头文件，编译器无法识别函数原型，引发未定义行为。
3. 多余的引用计数操作：Py_INCREF(inputArrayObj)是不必要的，输入参数的引用计数由Python解释器管理，额外增加会导致内存泄漏。

修复后的代码 #

#include <Python.h>
#include <stdio.h>
#include <stdlib.h>
#include <unistd.h>
#include <math.h>  // 新增：包含math.h头文件
#include <numpy/arrayobject.h>


// 修改参数类型为double*，内部变量也改为double
double calculateStandardDeviation(double *arr, int startIndex, int windowSize) {
    double mean = 0.0;
    double variance = 0.0;

    // Calculate the mean of the subarray
    for (int i = startIndex; i < startIndex + windowSize; i++) {
        mean += arr[i];
    }
    mean /= windowSize;

    // Calculate the variance
    for (int i = startIndex; i < startIndex + windowSize; i++) {
        double diff = arr[i] - mean;
        variance += diff * diff;  // 用乘法替代pow，更高效且避免函数调用开销
    }
    variance /= windowSize;

    // Calculate the standard deviation
    double standardDeviation = sqrt(variance);
    return standardDeviation;
}

static PyArrayObject *method_stdarray(PyObject * self, PyObject *args){
    PyArrayObject *inputArrayObj;
    int windowSize, stepSize;

    if (!PyArg_ParseTuple(args, "Oii", &inputArrayObj, &windowSize, &stepSize)) {
        return NULL;
    }
    if (!PyArray_Check(inputArrayObj)) {
        PyErr_SetString(PyExc_TypeError, "Input array must be numpy array.");
        return NULL;
    }
    int ndim = PyArray_NDIM(inputArrayObj);
    if (ndim > 1){
        PyErr_SetString(PyExc_BaseException, "Input numpy array has too many dimensions."); 
        return NULL;
    }
    if (PyArray_TYPE(inputArrayObj) != NPY_DOUBLE) {
        PyErr_SetString(PyExc_ValueError, "Input array must be of type float64.");
        return NULL;
    }

    int arraySize = PyArray_Size(inputArrayObj);
    double* inputArray = PyArray_DATA(inputArrayObj);

    if (inputArray == NULL) {
        PyErr_NoMemory();
        return NULL;
    }

    int numWindows = (arraySize - windowSize) / stepSize + 1;
    npy_intp numpyWindows = numWindows;

    PyArrayObject* resultArray = (PyArrayObject*)PyArray_EMPTY(1, &numpyWindows, NPY_DOUBLE, 0);
    if (resultArray == NULL) {
        PyErr_SetString(PyExc_MemoryError, "Failed to allocate memory for the result array.");
        return NULL;
    }

    double* resultData = PyArray_DATA(resultArray);
    
    for (int i = 0; i < numWindows; i++) {
        int startIndex = i * stepSize;
        // 修改为double类型接收结果
        double result = calculateStandardDeviation(inputArray, startIndex, windowSize);
        resultData[i] = result;
    }

    // 移除多余的Py_INCREF(inputArrayObj)

    return (PyObject*)resultArray;
}


static PyMethodDef meanyMethods[] = {
    {"vararray", method_stdarray, METH_VARARGS, "Takes a numpy array and two integers, returns a new array of standard deviations of subwindows"},
    {NULL, NULL, 0, NULL}
};


static struct PyModuleDef meanymodule = {
    PyModuleDef_HEAD_INIT,
    "meany",
    "My set of functions for analysing mea data",
    -1,
    meanyMethods
};

PyMODINIT_FUNC PyInit_meany(void) {
    import_array();
    return PyModule_Create(&meanymodule);
}

额外优化说明 #

• 将pow(arr[i]-mean, 2)替换为(arr[i]-mean)*(arr[i]-mean)，避免调用数学函数的开销，同时减少潜在的精度问题。
• 更新了方法注释，使其更准确描述功能。

内容的提问来源于stack exchange，提问作者FilterFeeder