139 |
gsl_vector_view xty; |
gsl_vector_view xty; |
140 |
gsl_vector_view xi; |
gsl_vector_view xi; |
141 |
gsl_vector_view xj; |
gsl_vector_view xj; |
142 |
|
gsl_vector *param_estimates; |
143 |
|
|
144 |
size_t i; |
size_t i; |
145 |
size_t j; |
size_t j; |
301 |
/* |
/* |
302 |
Use QR decomposition via GSL. |
Use QR decomposition via GSL. |
303 |
*/ |
*/ |
304 |
|
|
305 |
|
param_estimates = gsl_vector_alloc (1 + X->size2); |
306 |
design = gsl_matrix_alloc (X->size1, 1 + X->size2); |
design = gsl_matrix_alloc (X->size1, 1 + X->size2); |
307 |
|
|
308 |
for (j = 0; j < X->size1; j++) |
for (j = 0; j < X->size1; j++) |
316 |
} |
} |
317 |
gsl_multifit_linear_workspace *wk = |
gsl_multifit_linear_workspace *wk = |
318 |
gsl_multifit_linear_alloc (design->size1, design->size2); |
gsl_multifit_linear_alloc (design->size1, design->size2); |
319 |
rc = gsl_multifit_linear (design, Y, cache->param_estimates, |
rc = gsl_multifit_linear (design, Y, param_estimates, |
320 |
cache->cov, &(cache->sse), wk); |
cache->cov, &(cache->sse), wk); |
321 |
|
for (i = 0; i < cache->n_coeffs; i++) |
322 |
|
{ |
323 |
|
cache->coeff[i].estimate = gsl_vector_get (param_estimates, i); |
324 |
|
} |
325 |
if (rc == GSL_SUCCESS) |
if (rc == GSL_SUCCESS) |
326 |
{ |
{ |
327 |
gsl_multifit_linear_free (wk); |
gsl_multifit_linear_free (wk); |
328 |
|
gsl_vector_free (param_estimates); |
329 |
} |
} |
330 |
else |
else |
331 |
{ |
{ |