< Minuit | Indice | Appendici >
Strategia da seguire ed errori più comuni
Il modo migliore per imparare ad usare fminuit come fit non lineare consite nell'eseguire un esercizio.
L'esercizio che proponiamo richiede di generare dati con una legge nota ed aggiungere rumore normale, simulando così un esperimento. Successivamente si daranno i dati in pasto a fminuit per cercare di riottenerne i valori dei parametri della legge scelta. Illustriamo più in dettaglio i passi logici che occorre seguire nello svolgimento dell'esercizio:
- Scegliere una opportuna legge sperimentale, ad esempio:
{$ (1) \qquad\qquad y = A sin(\omega t + \varphi) e^{-t/T} $}
nella quale {$t$} è una variabile indipendente e {$y$} il valore misurato per un certo valore della variabile indipendente
- Generare dati sperimentali secondo questa legge, ossia scegliere un insieme discreto di valori di {$t$}, {$\{t_i\}$}, ad esempio equispaziati, per i quali si ottengono valori corrispondenti {$\{y_i\}$}, secondo la legge (1) e senza rumore. Aggiungere a questi dati un rumore normale di deviazione standard prefissata {$\sigma$}. Per generare questi dati occorre avere fissato il valore di cinque costanti, {$A, \omega, \phi, T$} che sono le costanti della legge, e {$\sigma$} che è l'errore delle misure.
- Attenzione! Controllare, visualizzando i dati ottenuti, di aver generato una legge sensata. Ad esempio, non è sensato usare una frequenza troppo alta rispetto al campionamento scelto. Ossia, se l'intervallo tra un tempo e il successivo è 1 s, non ha senso utilizzare frequenze maggiori di 0.5 Hz e conviene scegliere valori ben al di sotto di questo limite. Infatti la frequenza di 0.5 Hz implica un periodo di 2 s, corrispondenti a due soli punti sulla curva (il minimo assoluto per vedere direttamente l'oscillazione). All'estremo opposto occorre che la frequenza non abbia un periodo molto più lungo della lunghezza totale dell'intervallo scelto. Per calcolare rapidamente questi limiti opposti e scegliere valori intermedi conviene utilizzare per la funzione seno un argomento che esplicita la frequenza:
{$ (2) \qquad\qquad 2 \pi(\nu t + \frac {\varphi} {360}) $}
con {$\phi$} in gradi e {$\nu$} in Hz.
- Costruire una funzione {$\chi^2$} per la minimizzazione con MINUIT. Questa funzione può avere la seguente struttura:
function f=legge(p,d)
% commenti esplicativi
t=d(1,;);
A=p(1);w=p(2);phi=p(3);T=p(4);
f=A*sin(w*t+phi).*exp(-t/T);
if (size(d,1)==2), %chi-square, error = 1
f = sum((d(2,:) - f).^2);
elseif (size(d,1)>2), %chi-square, error = 3rd row of data
f = sum(((d(2,:) - f)./d(3,:)).^2);
end
- Prima di lanciare il fit, controllare di aver effettivamente generato i dati con la stessa legge con cui lo si esegue (ossia la
function legge): ciò si ottiene facendo disegnare dati e funzione alla funzione ausiliaria mnplot, fornita con fminuit, che consente di utilizzare proprio la stessa function legge per produrre il grafico. mnplot richiede in input i dati (una matrice con {$t$} nella prima riga e {$y$} nella seconda), i parametri (un vettore), e nome della funzione ('legge'). Se si utilizzano i parametri con cui si sono generati i dati i grafici dei punti e la curva della funzione devono coincidere entro l'errore {$\sigma$} aggiunto stocasticamente ai primi. Ecco un esempio di comadi:
dati=[t;y];
par=[A omega phi T];
mnplot(dati,par,'legge');
- Lo scopo del fit è di determinare ivalori delle quattro costanti {$A, \omega, \phi, T$} che meglio si adattano ai dati, e le relative incertezze.
fminuit e mnplot devono essere nel path di matlab. Per sincerarsene dare il comando help fminuit: deve comparire una descrizione del comando in inglese. Provare anche con help mnplot. Se fallisce solo il secondo vuol dire che si è aggiunto il solo path di fminuit senza le sue sottodirectories, che viceversa vanno incluse. Vedere istruzioni per l'istallazione
- A Minuit vanno passati il nome della funzione, il nome dalla funzione grafica d'appoggio, i parametri iniziali (diversi dai valori di
par, che in un vero esperimento non sarebbero noti, ma non troppo diversi, per evitare che il fit non li riesca a trovare) ed infine i dati (che devono includere una terza riga di errori). Ad esempio:
dati=[dati; sigma*ones(size(t))];
par=par*1.4;
[p,ep,c]=fminuit('legge','mnplot',par0,dati)
- Se non sono stati commessi errori comparirà il prompt di MINUIT, ossia la scritta
ENTER MINUIT COMMAND:
- Una sequenza standard di comandi sarà:
CALL 5. (che produce il grafico della funzione nel punto iniziale par0)
MIGR (che cerca il minimo col metodo di Newton)
CALL 5. (per controllare che il minimo sia stato trovato)
EXIT (per terminare)
- Un errore comune è quello di separare {$x, y$} ed errori, nella matrice
dati, con una virgola invece che con il punto e virgola, producendo un unico vettore. In questo caso la chiamata a mnplot risponde Nothing to plot
- Per ulteriori dettagli su MINUIT si rimanda al suo manuale
< Minuit | Indice | Appendici >
