tessl/pypi-python-casacore

Python bindings for the CASACORE radio astronomy library providing comprehensive interfaces for table operations, image processing, coordinate systems, and astronomical measurements.

—

Pending

Overview

Eval results

Files

Fitting Operations

Name: tessl/pypi-python-casacore
Author: tessl

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.

Core Imports

from casacore.fitting import fitserver

Capabilities

Fit Server

Central interface for least squares fitting with support for multiple simultaneous fits and constraint handling.

class fitserver:
    def __init__(self, **kwargs):
        """
        Create fitting server.
        
        Parameters:
        - ftype: str, fitting type ('real' or 'complex', default 'real')
        - constraint_type: str, constraint handling method
        """
    
    def fitter(self, **kwargs):
        """
        Create sub-fitter instance.
        
        Parameters:
        - Same as __init__
        
        Returns:
        int, fitter ID for use with fid parameter
        """
    
    def set(self, **kwargs):
        """
        Set fitting parameters.
        
        Parameters:
        - ftype: str, fitting type ('real' or 'complex')
        - constraint_type: str, constraint method
        - fid: int, fitter ID (optional)
        """

Linear Fitting

Solve linear least squares problems with optional constraints and weights.

class fitserver:
    def linear(self, func, x, y, **kwargs):
        """
        Perform linear least squares fit.
        
        Parameters:
        - func: functional object, model function (must be linear in parameters)
        - x: numpy array, independent variable(s) 
        - y: numpy array, dependent variable (observations)
        - weight: numpy array, weights for each observation (optional)
        - sigma: numpy array, standard deviations (optional, alternative to weight)
        - fid: int, fitter ID for sub-fitter (optional)
        
        The functional must be linear in its parameters for this method.
        For non-linear problems, use functional() method.
        
        Examples of linear functionals:
        - polynomial: p0 + p1*x + p2*x²
        - sinusoid: p0*sin(x) + p1*cos(x)  
        - exponential: p0*exp(x) + p1*exp(2*x)
        
        Examples of non-linear functionals:
        - Gaussian: p0*exp(-((x-p1)/p2)²)  [non-linear in p1, p2]
        - exponential: p0*exp(p1*x)        [non-linear in p1]
        """

Non-Linear Fitting

Solve non-linear least squares problems using iterative algorithms.

class fitserver:
    def functional(self, func, x, y, **kwargs):
        """
        Perform non-linear least squares fit.
        
        Parameters:
        - func: functional object, model function (can be non-linear)
        - x: numpy array, independent variable(s)
        - y: numpy array, dependent variable (observations)
        - weight: numpy array, observation weights (optional)
        - sigma: numpy array, standard deviations (optional)
        - maxiter: int, maximum iterations (default 1000)
        - criteria: float, convergence criterion (default 0.001)
        - fid: int, fitter ID (optional)
        
        The functional can be non-linear in parameters. Initial parameter
        values should be set in the functional before calling this method.
        Good initial guesses are important for convergence.
        
        Algorithm: Levenberg-Marquardt with automatic derivative calculation
        """

Solution Analysis

Extract fitting results including parameter values, errors, and quality metrics.

class fitserver:
    def solution(self, fid=None):
        """
        Get fitted parameter values.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        numpy array, optimized parameter values
        """
    
    def error(self, fid=None):
        """
        Get parameter errors (standard deviations).
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        numpy array, parameter standard deviations
        """
    
    def covariance(self, fid=None):
        """
        Get parameter covariance matrix.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        2D numpy array, covariance matrix
        """
    
    def correlation(self, fid=None):
        """
        Get parameter correlation matrix.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        2D numpy array, correlation coefficients (-1 to 1)
        """
    
    def chi2(self, fid=None):
        """
        Get chi-squared statistic.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        float, chi-squared value
        """
    
    def sd(self, fid=None):
        """
        Get standard deviation per observation.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        float, RMS residual
        """
    
    def stddev(self, fid=None):
        """
        Get standard deviation per unit weight.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        float, weighted RMS
        """
    
    def rank(self, fid=None):
        """
        Get solution rank (number of independent parameters).
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        int, matrix rank
        """
    
    def deficiency(self, fid=None):
        """
        Get rank deficiency (nparams - rank).
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        int, rank deficiency
        """

Constraint Management

Add and manage linear constraints on fitted parameters.

class fitserver:
    def addconstraint(self, x=[], y=0.0, **kwargs):
        """
        Add linear constraint to fit.
        
        Parameters:
        - x: list or numpy array, constraint coefficients
        - y: float, constraint value (default 0.0)
        - fid: int, fitter ID (optional)
        
        Constraint form: sum(x[i] * param[i]) = y
        
        Examples:
        # Sum of parameters equals 1.0
        addconstraint([1, 1, 1], 1.0)
        
        # Second parameter fixed at 5.0  
        addconstraint([0, 1, 0], 5.0)
        
        # Difference between param[0] and param[1] equals 0
        addconstraint([1, -1], 0.0)
        """
    
    def clearconstraints(self, fid=None):
        """
        Remove all constraints.
        
        Parameters:
        - fid: int, fitter ID (optional)
        """
    
    def nconstraints(self, fid=None):
        """
        Get number of active constraints.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        int, number of constraints
        """

SVD and Rank-Deficient Solutions

Handle rank-deficient systems using Singular Value Decomposition.

class fitserver:
    def constraint(self, which=0, fid=None):
        """
        Get SVD constraint for rank-deficient system.
        
        Parameters:
        - which: int, constraint index (0-based)
        - fid: int, fitter ID (optional)
        
        Returns:
        numpy array, constraint coefficients
        
        When a fit is rank-deficient, SVD provides constraint equations
        that can be used to obtain a unique solution.
        """
    
    def svd(self, fid=None):
        """
        Get SVD decomposition information.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        dict with SVD matrices and singular values
        """

Complex Fitting

Handle complex-valued data and parameters with specialized algorithms.

class fitserver:
    def set_complex(self, fid=None):
        """
        Set fitter to complex mode.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        In complex mode:
        - Data can be complex-valued
        - Parameters can be complex
        - Constraints can involve complex relationships
        - Conjugate constraints are supported
        """
    
    def real(self, fid=None):
        """
        Set fitter to real mode (default).
        
        Parameters:
        - fid: int, fitter ID (optional)
        """
    
    def complex(self, fid=None):
        """
        Set fitter to complex mode.
        
        Parameters:  
        - fid: int, fitter ID (optional)
        """

Fit Quality Assessment

Evaluate fit quality and detect potential problems.

class fitserver:
    def residuals(self, fid=None):
        """
        Get fit residuals (observed - model).
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        numpy array, residual values
        """
    
    def fitted(self, fid=None):
        """
        Get fitted model values.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        numpy array, model predictions at data points
        """
    
    def reducedchi2(self, fid=None):
        """
        Get reduced chi-squared (chi²/dof).
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        float, reduced chi-squared statistic
        """
    
    def dof(self, fid=None):
        """
        Get degrees of freedom.
        
        Parameters:
        - fid: int, fitter ID (optional)
        
        Returns:
        int, degrees of freedom (ndata - nparams + nconstraints)
        """

Usage Examples

Basic Linear Fitting

from casacore.fitting import fitserver
from casacore.functionals import poly, gaussian1d
import numpy as np

# Create synthetic data
x = np.linspace(0, 10, 50)
y_true = 2.0 + 1.5*x + 0.3*x**2  # True quadratic
noise = np.random.normal(0, 0.5, len(x))
y = y_true + noise

# Create fitting server and polynomial model
fit = fitserver()
p = poly(2)  # 2nd order polynomial: p0 + p1*x + p2*x²

# Perform linear fit
fit.linear(p, x, y)

# Get results
params = fit.solution()
errors = fit.error()
chi2 = fit.chi2()
rms = fit.sd()

print(f"Fitted parameters: {params}")
print(f"Parameter errors: {errors}")
print(f"Chi-squared: {chi2:.3f}")
print(f"RMS residual: {rms:.3f}")
print(f"True parameters: [2.0, 1.5, 0.3]")

Non-Linear Fitting

# Create synthetic Gaussian data
x = np.linspace(-5, 15, 100)
y_true = 10.0 * np.exp(-0.5*((x - 5.0)/2.0)**2) + 1.0
noise = np.random.normal(0, 0.2, len(x))
y = y_true + noise

# Create Gaussian functional with initial guess
gauss = gaussian1d(peak=8.0, center=4.0, width=1.5)  # Initial guess

# Perform non-linear fit
fit.functional(gauss, x, y)

# Get results
params = fit.solution()
errors = fit.error()
covar = fit.covariance()

print(f"Fitted Gaussian:")
print(f"  Peak: {params[0]:.3f} ± {errors[0]:.3f}")
print(f"  Center: {params[1]:.3f} ± {errors[1]:.3f}")  
print(f"  Width: {params[2]:.3f} ± {errors[2]:.3f}")
print(f"True parameters: [10.0, 5.0, 2.0]")

# Calculate fitted curve
y_fit = gauss(x)
plt.plot(x, y, 'o', label='Data')
plt.plot(x, y_fit, '-', label='Fit')
plt.legend()

Constrained Fitting

# Fit triangle angles with constraint that sum = 180°
angles_data = np.array([49.5, 60.8, 69.2])  # Measured angles
angle_errors = np.array([0.5, 0.5, 0.5])    # Measurement errors

# Create linear model for three independent parameters
from casacore.functionals import compiled
model = compiled('p*x + p1*x1 + p2*x2')  # p0*I1 + p1*I2 + p2*I3

# Design matrix (identity for this simple case)
x_design = np.array([[1, 0, 0],  # First angle
                     [0, 1, 0],  # Second angle  
                     [0, 0, 1]]) # Third angle

# Fit without constraints
fit.linear(model, x_design.T, angles_data, sigma=angle_errors)
unconstrained = fit.solution()
print(f"Unconstrained fit: {unconstrained}, sum = {np.sum(unconstrained):.1f}°")

# Add constraint: sum of angles = 180°
fit.addconstraint([1, 1, 1], 180.0)
fit.linear(model, x_design.T, angles_data, sigma=angle_errors)
constrained = fit.solution()
print(f"Constrained fit: {constrained}, sum = {np.sum(constrained):.1f}°")

# Errors are now correlated due to constraint
errors_constrained = fit.error()
correlation = fit.correlation()
print(f"Parameter errors: {errors_constrained}")
print(f"Correlation matrix:\n{correlation}")

Multiple Fits with Sub-Fitters

# Fit multiple datasets simultaneously
datasets = [
    (x1, y1),  # First dataset
    (x2, y2),  # Second dataset  
    (x3, y3)   # Third dataset
]

# Create sub-fitters
fitter_ids = []
for i in range(len(datasets)):
    fid = fit.fitter()  # Create new sub-fitter
    fitter_ids.append(fid)

# Fit each dataset
results = []
for i, (x_data, y_data) in enumerate(datasets):
    gauss = gaussian1d(peak=5.0, center=0.0, width=1.0)  # Initial guess
    fit.functional(gauss, x_data, y_data, fid=fitter_ids[i])
    
    params = fit.solution(fid=fitter_ids[i])
    errors = fit.error(fid=fitter_ids[i])
    chi2 = fit.chi2(fid=fitter_ids[i])
    
    results.append({
        'params': params,
        'errors': errors, 
        'chi2': chi2
    })

for i, result in enumerate(results):
    print(f"Dataset {i+1}: peak={result['params'][0]:.2f}, χ²={result['chi2']:.2f}")

Complex Data Fitting

# Fit complex-valued data (e.g., visibility amplitudes and phases)
freq = np.linspace(1.4, 1.5, 50)  # GHz
# Complex visibility with amplitude and phase structure
visibility = (5.0 + 2.0j) * np.exp(-((freq - 1.42)/0.05)**2) * np.exp(1j * 0.1 * freq)
noise_real = np.random.normal(0, 0.1, len(freq))
noise_imag = np.random.normal(0, 0.1, len(freq))
vis_noisy = visibility + noise_real + 1j * noise_imag

# Create complex fitter
complex_fit = fitserver(ftype='complex')

# Create complex Gaussian model  
from casacore.functionals import compiled
complex_model = compiled('(p + p1*1j) * exp(-((x - p2)/p3)^2) * exp(1j*p4*x)')
# Parameters: [real_amp, imag_amp, center_freq, width, phase_slope]

# Set initial parameters
complex_model.set_parameters([4.0, 1.5, 1.42, 0.05, 0.1])

# Fit complex data
complex_fit.functional(complex_model, freq, vis_noisy)

# Get complex results
params = complex_fit.solution()
errors = complex_fit.error()

print(f"Complex amplitude: {params[0]:.2f} + {params[1]:.2f}j")
print(f"Center frequency: {params[2]:.4f} GHz")
print(f"Spectral width: {params[3]:.4f} GHz")
print(f"Phase slope: {params[4]:.3f} rad/GHz")

Rank-Deficient Problems and SVD

# Create rank-deficient system (more parameters than constraints)
x = np.linspace(0, 1, 10)
y = np.ones(10)  # Constant data

# Try to fit with too many parameters
high_order_poly = poly(5)  # 6 parameters for constant data
fit.linear(high_order_poly, x, y)

# Check rank deficiency
rank = fit.rank()
deficiency = fit.deficiency()
nparams = high_order_poly.nparameters()

print(f"Parameters: {nparams}, Rank: {rank}, Deficiency: {deficiency}")

if deficiency > 0:
    print("System is rank-deficient. SVD constraints available:")
    for i in range(deficiency):
        constraint = fit.constraint(i)
        print(f"  Constraint {i}: {constraint}")
    
    # Add SVD constraints to get unique solution
    for i in range(deficiency):
        svd_constraint = fit.constraint(i)
        fit.addconstraint(svd_constraint, 0.0)
    
    # Refit with constraints
    fit.linear(high_order_poly, x, y)
    solution = fit.solution()
    print(f"Constrained solution: {solution}")

Install with Tessl CLI