0
# Fitting Operations
1
2
Linear and non-linear least squares fitting with constraint support. Handles both real and complex fitting problems with SVD solutions for rank-deficient systems, comprehensive error analysis, and advanced constraint mechanisms for scientific data analysis.
3
4
## Core Imports
5
6
```python
7
from casacore.fitting import fitserver
8
```
9
10
## Capabilities
11
12
### Fit Server
13
14
Central interface for least squares fitting with support for multiple simultaneous fits and constraint handling.
15
16
```python { .api }
17
class fitserver:
18
def __init__(self, **kwargs):
19
"""
20
Create fitting server.
21
22
Parameters:
23
- ftype: str, fitting type ('real' or 'complex', default 'real')
24
- constraint_type: str, constraint handling method
25
"""
26
27
def fitter(self, **kwargs):
28
"""
29
Create sub-fitter instance.
30
31
Parameters:
32
- Same as __init__
33
34
Returns:
35
int, fitter ID for use with fid parameter
36
"""
37
38
def set(self, **kwargs):
39
"""
40
Set fitting parameters.
41
42
Parameters:
43
- ftype: str, fitting type ('real' or 'complex')
44
- constraint_type: str, constraint method
45
- fid: int, fitter ID (optional)
46
"""
47
```
48
49
### Linear Fitting
50
51
Solve linear least squares problems with optional constraints and weights.
52
53
```python { .api }
54
class fitserver:
55
def linear(self, func, x, y, **kwargs):
56
"""
57
Perform linear least squares fit.
58
59
Parameters:
60
- func: functional object, model function (must be linear in parameters)
61
- x: numpy array, independent variable(s)
62
- y: numpy array, dependent variable (observations)
63
- weight: numpy array, weights for each observation (optional)
64
- sigma: numpy array, standard deviations (optional, alternative to weight)
65
- fid: int, fitter ID for sub-fitter (optional)
66
67
The functional must be linear in its parameters for this method.
68
For non-linear problems, use functional() method.
69
70
Examples of linear functionals:
71
- polynomial: p0 + p1*x + p2*x²
72
- sinusoid: p0*sin(x) + p1*cos(x)
73
- exponential: p0*exp(x) + p1*exp(2*x)
74
75
Examples of non-linear functionals:
76
- Gaussian: p0*exp(-((x-p1)/p2)²) [non-linear in p1, p2]
77
- exponential: p0*exp(p1*x) [non-linear in p1]
78
"""
79
```
80
81
### Non-Linear Fitting
82
83
Solve non-linear least squares problems using iterative algorithms.
84
85
```python { .api }
86
class fitserver:
87
def functional(self, func, x, y, **kwargs):
88
"""
89
Perform non-linear least squares fit.
90
91
Parameters:
92
- func: functional object, model function (can be non-linear)
93
- x: numpy array, independent variable(s)
94
- y: numpy array, dependent variable (observations)
95
- weight: numpy array, observation weights (optional)
96
- sigma: numpy array, standard deviations (optional)
97
- maxiter: int, maximum iterations (default 1000)
98
- criteria: float, convergence criterion (default 0.001)
99
- fid: int, fitter ID (optional)
100
101
The functional can be non-linear in parameters. Initial parameter
102
values should be set in the functional before calling this method.
103
Good initial guesses are important for convergence.
104
105
Algorithm: Levenberg-Marquardt with automatic derivative calculation
106
"""
107
```
108
109
### Solution Analysis
110
111
Extract fitting results including parameter values, errors, and quality metrics.
112
113
```python { .api }
114
class fitserver:
115
def solution(self, fid=None):
116
"""
117
Get fitted parameter values.
118
119
Parameters:
120
- fid: int, fitter ID (optional)
121
122
Returns:
123
numpy array, optimized parameter values
124
"""
125
126
def error(self, fid=None):
127
"""
128
Get parameter errors (standard deviations).
129
130
Parameters:
131
- fid: int, fitter ID (optional)
132
133
Returns:
134
numpy array, parameter standard deviations
135
"""
136
137
def covariance(self, fid=None):
138
"""
139
Get parameter covariance matrix.
140
141
Parameters:
142
- fid: int, fitter ID (optional)
143
144
Returns:
145
2D numpy array, covariance matrix
146
"""
147
148
def correlation(self, fid=None):
149
"""
150
Get parameter correlation matrix.
151
152
Parameters:
153
- fid: int, fitter ID (optional)
154
155
Returns:
156
2D numpy array, correlation coefficients (-1 to 1)
157
"""
158
159
def chi2(self, fid=None):
160
"""
161
Get chi-squared statistic.
162
163
Parameters:
164
- fid: int, fitter ID (optional)
165
166
Returns:
167
float, chi-squared value
168
"""
169
170
def sd(self, fid=None):
171
"""
172
Get standard deviation per observation.
173
174
Parameters:
175
- fid: int, fitter ID (optional)
176
177
Returns:
178
float, RMS residual
179
"""
180
181
def stddev(self, fid=None):
182
"""
183
Get standard deviation per unit weight.
184
185
Parameters:
186
- fid: int, fitter ID (optional)
187
188
Returns:
189
float, weighted RMS
190
"""
191
192
def rank(self, fid=None):
193
"""
194
Get solution rank (number of independent parameters).
195
196
Parameters:
197
- fid: int, fitter ID (optional)
198
199
Returns:
200
int, matrix rank
201
"""
202
203
def deficiency(self, fid=None):
204
"""
205
Get rank deficiency (nparams - rank).
206
207
Parameters:
208
- fid: int, fitter ID (optional)
209
210
Returns:
211
int, rank deficiency
212
"""
213
```
214
215
### Constraint Management
216
217
Add and manage linear constraints on fitted parameters.
218
219
```python { .api }
220
class fitserver:
221
def addconstraint(self, x=[], y=0.0, **kwargs):
222
"""
223
Add linear constraint to fit.
224
225
Parameters:
226
- x: list or numpy array, constraint coefficients
227
- y: float, constraint value (default 0.0)
228
- fid: int, fitter ID (optional)
229
230
Constraint form: sum(x[i] * param[i]) = y
231
232
Examples:
233
# Sum of parameters equals 1.0
234
addconstraint([1, 1, 1], 1.0)
235
236
# Second parameter fixed at 5.0
237
addconstraint([0, 1, 0], 5.0)
238
239
# Difference between param[0] and param[1] equals 0
240
addconstraint([1, -1], 0.0)
241
"""
242
243
def clearconstraints(self, fid=None):
244
"""
245
Remove all constraints.
246
247
Parameters:
248
- fid: int, fitter ID (optional)
249
"""
250
251
def nconstraints(self, fid=None):
252
"""
253
Get number of active constraints.
254
255
Parameters:
256
- fid: int, fitter ID (optional)
257
258
Returns:
259
int, number of constraints
260
"""
261
```
262
263
### SVD and Rank-Deficient Solutions
264
265
Handle rank-deficient systems using Singular Value Decomposition.
266
267
```python { .api }
268
class fitserver:
269
def constraint(self, which=0, fid=None):
270
"""
271
Get SVD constraint for rank-deficient system.
272
273
Parameters:
274
- which: int, constraint index (0-based)
275
- fid: int, fitter ID (optional)
276
277
Returns:
278
numpy array, constraint coefficients
279
280
When a fit is rank-deficient, SVD provides constraint equations
281
that can be used to obtain a unique solution.
282
"""
283
284
def svd(self, fid=None):
285
"""
286
Get SVD decomposition information.
287
288
Parameters:
289
- fid: int, fitter ID (optional)
290
291
Returns:
292
dict with SVD matrices and singular values
293
"""
294
```
295
296
### Complex Fitting
297
298
Handle complex-valued data and parameters with specialized algorithms.
299
300
```python { .api }
301
class fitserver:
302
def set_complex(self, fid=None):
303
"""
304
Set fitter to complex mode.
305
306
Parameters:
307
- fid: int, fitter ID (optional)
308
309
In complex mode:
310
- Data can be complex-valued
311
- Parameters can be complex
312
- Constraints can involve complex relationships
313
- Conjugate constraints are supported
314
"""
315
316
def real(self, fid=None):
317
"""
318
Set fitter to real mode (default).
319
320
Parameters:
321
- fid: int, fitter ID (optional)
322
"""
323
324
def complex(self, fid=None):
325
"""
326
Set fitter to complex mode.
327
328
Parameters:
329
- fid: int, fitter ID (optional)
330
"""
331
```
332
333
### Fit Quality Assessment
334
335
Evaluate fit quality and detect potential problems.
336
337
```python { .api }
338
class fitserver:
339
def residuals(self, fid=None):
340
"""
341
Get fit residuals (observed - model).
342
343
Parameters:
344
- fid: int, fitter ID (optional)
345
346
Returns:
347
numpy array, residual values
348
"""
349
350
def fitted(self, fid=None):
351
"""
352
Get fitted model values.
353
354
Parameters:
355
- fid: int, fitter ID (optional)
356
357
Returns:
358
numpy array, model predictions at data points
359
"""
360
361
def reducedchi2(self, fid=None):
362
"""
363
Get reduced chi-squared (chi²/dof).
364
365
Parameters:
366
- fid: int, fitter ID (optional)
367
368
Returns:
369
float, reduced chi-squared statistic
370
"""
371
372
def dof(self, fid=None):
373
"""
374
Get degrees of freedom.
375
376
Parameters:
377
- fid: int, fitter ID (optional)
378
379
Returns:
380
int, degrees of freedom (ndata - nparams + nconstraints)
381
"""
382
```
383
384
## Usage Examples
385
386
### Basic Linear Fitting
387
388
```python
389
from casacore.fitting import fitserver
390
from casacore.functionals import poly, gaussian1d
391
import numpy as np
392
393
# Create synthetic data
394
x = np.linspace(0, 10, 50)
395
y_true = 2.0 + 1.5*x + 0.3*x**2 # True quadratic
396
noise = np.random.normal(0, 0.5, len(x))
397
y = y_true + noise
398
399
# Create fitting server and polynomial model
400
fit = fitserver()
401
p = poly(2) # 2nd order polynomial: p0 + p1*x + p2*x²
402
403
# Perform linear fit
404
fit.linear(p, x, y)
405
406
# Get results
407
params = fit.solution()
408
errors = fit.error()
409
chi2 = fit.chi2()
410
rms = fit.sd()
411
412
print(f"Fitted parameters: {params}")
413
print(f"Parameter errors: {errors}")
414
print(f"Chi-squared: {chi2:.3f}")
415
print(f"RMS residual: {rms:.3f}")
416
print(f"True parameters: [2.0, 1.5, 0.3]")
417
```
418
419
### Non-Linear Fitting
420
421
```python
422
# Create synthetic Gaussian data
423
x = np.linspace(-5, 15, 100)
424
y_true = 10.0 * np.exp(-0.5*((x - 5.0)/2.0)**2) + 1.0
425
noise = np.random.normal(0, 0.2, len(x))
426
y = y_true + noise
427
428
# Create Gaussian functional with initial guess
429
gauss = gaussian1d(peak=8.0, center=4.0, width=1.5) # Initial guess
430
431
# Perform non-linear fit
432
fit.functional(gauss, x, y)
433
434
# Get results
435
params = fit.solution()
436
errors = fit.error()
437
covar = fit.covariance()
438
439
print(f"Fitted Gaussian:")
440
print(f" Peak: {params[0]:.3f} ± {errors[0]:.3f}")
441
print(f" Center: {params[1]:.3f} ± {errors[1]:.3f}")
442
print(f" Width: {params[2]:.3f} ± {errors[2]:.3f}")
443
print(f"True parameters: [10.0, 5.0, 2.0]")
444
445
# Calculate fitted curve
446
y_fit = gauss(x)
447
plt.plot(x, y, 'o', label='Data')
448
plt.plot(x, y_fit, '-', label='Fit')
449
plt.legend()
450
```
451
452
### Constrained Fitting
453
454
```python
455
# Fit triangle angles with constraint that sum = 180°
456
angles_data = np.array([49.5, 60.8, 69.2]) # Measured angles
457
angle_errors = np.array([0.5, 0.5, 0.5]) # Measurement errors
458
459
# Create linear model for three independent parameters
460
from casacore.functionals import compiled
461
model = compiled('p*x + p1*x1 + p2*x2') # p0*I1 + p1*I2 + p2*I3
462
463
# Design matrix (identity for this simple case)
464
x_design = np.array([[1, 0, 0], # First angle
465
[0, 1, 0], # Second angle
466
[0, 0, 1]]) # Third angle
467
468
# Fit without constraints
469
fit.linear(model, x_design.T, angles_data, sigma=angle_errors)
470
unconstrained = fit.solution()
471
print(f"Unconstrained fit: {unconstrained}, sum = {np.sum(unconstrained):.1f}°")
472
473
# Add constraint: sum of angles = 180°
474
fit.addconstraint([1, 1, 1], 180.0)
475
fit.linear(model, x_design.T, angles_data, sigma=angle_errors)
476
constrained = fit.solution()
477
print(f"Constrained fit: {constrained}, sum = {np.sum(constrained):.1f}°")
478
479
# Errors are now correlated due to constraint
480
errors_constrained = fit.error()
481
correlation = fit.correlation()
482
print(f"Parameter errors: {errors_constrained}")
483
print(f"Correlation matrix:\n{correlation}")
484
```
485
486
### Multiple Fits with Sub-Fitters
487
488
```python
489
# Fit multiple datasets simultaneously
490
datasets = [
491
(x1, y1), # First dataset
492
(x2, y2), # Second dataset
493
(x3, y3) # Third dataset
494
]
495
496
# Create sub-fitters
497
fitter_ids = []
498
for i in range(len(datasets)):
499
fid = fit.fitter() # Create new sub-fitter
500
fitter_ids.append(fid)
501
502
# Fit each dataset
503
results = []
504
for i, (x_data, y_data) in enumerate(datasets):
505
gauss = gaussian1d(peak=5.0, center=0.0, width=1.0) # Initial guess
506
fit.functional(gauss, x_data, y_data, fid=fitter_ids[i])
507
508
params = fit.solution(fid=fitter_ids[i])
509
errors = fit.error(fid=fitter_ids[i])
510
chi2 = fit.chi2(fid=fitter_ids[i])
511
512
results.append({
513
'params': params,
514
'errors': errors,
515
'chi2': chi2
516
})
517
518
for i, result in enumerate(results):
519
print(f"Dataset {i+1}: peak={result['params'][0]:.2f}, χ²={result['chi2']:.2f}")
520
```
521
522
### Complex Data Fitting
523
524
```python
525
# Fit complex-valued data (e.g., visibility amplitudes and phases)
526
freq = np.linspace(1.4, 1.5, 50) # GHz
527
# Complex visibility with amplitude and phase structure
528
visibility = (5.0 + 2.0j) * np.exp(-((freq - 1.42)/0.05)**2) * np.exp(1j * 0.1 * freq)
529
noise_real = np.random.normal(0, 0.1, len(freq))
530
noise_imag = np.random.normal(0, 0.1, len(freq))
531
vis_noisy = visibility + noise_real + 1j * noise_imag
532
533
# Create complex fitter
534
complex_fit = fitserver(ftype='complex')
535
536
# Create complex Gaussian model
537
from casacore.functionals import compiled
538
complex_model = compiled('(p + p1*1j) * exp(-((x - p2)/p3)^2) * exp(1j*p4*x)')
539
# Parameters: [real_amp, imag_amp, center_freq, width, phase_slope]
540
541
# Set initial parameters
542
complex_model.set_parameters([4.0, 1.5, 1.42, 0.05, 0.1])
543
544
# Fit complex data
545
complex_fit.functional(complex_model, freq, vis_noisy)
546
547
# Get complex results
548
params = complex_fit.solution()
549
errors = complex_fit.error()
550
551
print(f"Complex amplitude: {params[0]:.2f} + {params[1]:.2f}j")
552
print(f"Center frequency: {params[2]:.4f} GHz")
553
print(f"Spectral width: {params[3]:.4f} GHz")
554
print(f"Phase slope: {params[4]:.3f} rad/GHz")
555
```
556
557
### Rank-Deficient Problems and SVD
558
559
```python
560
# Create rank-deficient system (more parameters than constraints)
561
x = np.linspace(0, 1, 10)
562
y = np.ones(10) # Constant data
563
564
# Try to fit with too many parameters
565
high_order_poly = poly(5) # 6 parameters for constant data
566
fit.linear(high_order_poly, x, y)
567
568
# Check rank deficiency
569
rank = fit.rank()
570
deficiency = fit.deficiency()
571
nparams = high_order_poly.nparameters()
572
573
print(f"Parameters: {nparams}, Rank: {rank}, Deficiency: {deficiency}")
574
575
if deficiency > 0:
576
print("System is rank-deficient. SVD constraints available:")
577
for i in range(deficiency):
578
constraint = fit.constraint(i)
579
print(f" Constraint {i}: {constraint}")
580
581
# Add SVD constraints to get unique solution
582
for i in range(deficiency):
583
svd_constraint = fit.constraint(i)
584
fit.addconstraint(svd_constraint, 0.0)
585
586
# Refit with constraints
587
fit.linear(high_order_poly, x, y)
588
solution = fit.solution()
589
print(f"Constrained solution: {solution}")
590
```