Я пытаюсь согласовать собранные данные с полиномиальным уравнением, и я нашел функцию lfit из числовых рецептов. У меня есть доступ только ко второму изданию, поэтому я использую это.
Я прочитал о функции lfit и ее параметрах, одним из которых является указатель функции, указанный в документации как
void (*funcs)(float, float [], int))
с помощью
Пользователь предоставляет подпрограмму funcs (x, afunc, ma), которая возвращает ma базовых функций, оцененных при x = x, в массиве afunc [1..ma].
Я пытаюсь понять, как работает эта lfit функция. Пример функции, которую я нашел, приведен ниже:
void fpoly(float x, float p[], int np)
/*Fitting routine for a polynomial of degree np-1, with coefficients in the array p[1..np].*/
{
int j;
p[1]=1.0;
for (j=2;j<=np;j++)
p[j]=p[j-1]*x;
}
Когда я запускаю исходный код функции lfit в gdb, я не вижу ссылки на указатель funcs. Когда я пытаюсь установить простой набор данных с помощью функции, я получаю следующее сообщение об ошибке.
Numerical Recipes run-time error...
gaussj: Singular Matrix
...now exiting to system...
Ясно, что матрица как-то определяется со всеми нулями.
Я собираюсь включить подгонку этой функции в большой цикл, поэтому использование другого языка на самом деле не вариант. Следовательно, почему я планирую использовать C / C ++.
Для справки, тестовая программа приведена здесь:
int main()
{
float x[5] = {0., 0., 1., 2., 3.};
float y[5] = {0., 0., 1.2, 3.9, 7.5};
float sig[5] = {1., 1., 1., 1., 1.};
int ndat = 4;
int ma = 4; /* parameters in equation */
float a[5] = {1, 1, 1, 0.1, 1.5};
int ia[5] = {1, 1, 1, 1, 1};
float **covar = matrix(1, ma, 1, ma);
float chisq = 0;
lfit(x,y,sig,ndat,a,ia,ma,covar,&chisq,fpoly);
printf("%f\n", chisq);
free_matrix(covar, 1, ma, 1, ma);
return 0;
}
Также запутывает проблему, все функции Числовых Рецептов индексируются 1 массивом, поэтому, если у кого-то есть исправления в моих объявлениях массива, сообщите мне также!
Приветствия
edit: Хорошо, только что обнаружил проблему с этим.
При копировании кода числовых рецептов я случайно закомментировал некоторый реальный код, думая, что это просто комментарий. Я больше не получаю единственное значение ошибки