AGGIORNAMENTO ESPERIMENTO 2:Nelle puntate precedenti...Nel primo esperimento avevo generato dei dati "sintetici" per mettere alla prova i vari algoritmi esistenti e alcuni nuovi algoritmi. Questi dati erano rappresentati da un segnale
gaussiano, tutti con lo stesso valore medio (125) e la stessa deviazione standard, a cui andavo ad aggiungere una certa percentuale (tra il 20 e il 40%) di dati
outliers, ovvero dati "anomali", rappresentati da una distribuzione
uniforme simmetrica rispetto alla media.
Allegato:
pdf.jpg [ 14.81 KiB | Osservato 1296 volte ]
I risultati erano stati interessanti, con il metodo sigma-clipping e l'autoadaptive weighted average che si comportavano sostanzialmente allo stesso modo, entrambi comunque non si comportavano particolarmente meglio della media.
Il metodo del fitting-gaussiano (quello su cui sto/stiamo lavorando) non dava risultati interessanti.
In questa puntataGrazie all'idea di Andrea Console ho avuto una idea: ho ipotizzato di avere un set di 50 immagini da 0.5MPixel (di più erano troppo pesanti da processare).
Innanzitutto ho migliorato la creazione dei dati sintetici. Come ha sottolineato Andrea, gli outliers non sono simmetrici. Ho quindi creato una distribuzione di outliers asimmetrica tra 125 e 255.
Allegato:
pdf2.jpg [ 13.28 KiB | Osservato 1296 volte ]
Per ognuno di questi pixel si può applicare uno dei classici algoritmi. Per ogni pixel, questi lavorano sui 50 valori che questo assume in ciascuna immagine. I risultati sono i soliti.
Qui però entra in gioco l'idea di Andrea. Poiché noi abbiamo 500'000 pixel x 50 foto, possiamo stimare abbastanza bene la distribuzione di probabilità dei dati di ciascun pixel a cui abbiamo sottratto preventivamente il suo valore medio.
A questo punto per ciascun pixel, per il quale abbiamo 50 dati, possiamo trovare il valore medio che minimizza l'errore rispetto alla distribuzione di probabilità stimata. In altre parole, ipotizziamo che per un numero di foto che tendesse all'infinito, il segnale dell'i-esimo pixel tenderebbe a quella che abbiamo stimato sull'intero sensore. Avendo però un numero finito di campioni, cerchiamo di stimare quale sia il valore medio a partire da una regressione ai minimi quadrati.
Affinchè questo metodo abbia senso occorre conoscere la risposta a due domande:
Domanda 1: nella realtà, se prendiamo due pixels che misurano due aree di cielo con la stessa intensità, la distribuzione di probabilità del segnale è uguale per tutti e due?
Domanda 2: nel caso precedente, se prendiamo due pixels che misurano due aree di cielo con intensità diversa, la distribuzione di probabilità del segnale è uguale per entrambe, una volta che si sottrae il valor medio?
La domanda 1 è una condizione necessaria affinché il metodo funzioni. Se la risposta è NO, possiamo anche fermarci qui.
La domanda 2 non è necessaria, ma aiuterebbe non poco. In caso di risposta negativa infatti non potremmo più limitarci a stimare una "master PDF", ma dovremmo stimarne una per ciascun valore di intensità media.
Detto questo, ho applicato il metodo ai dati sintetici di cui vi parlavo. Il metodo è molto lento per ora (non ho fatto nessuna ottimizzazione in questo senso) e ogni iterazione richiede circa 2 ore di calcolo. Per ora quindi ho fatto una sola prova.
I risultati sono incoraggianti.
Come potete vedere nell'immagine seguente la PDF dei dati viene stimata sufficiente bene già a partire dalla seconda iterazione del metodo. Le iterazioni successive migliorano leggermente la stima della moda del segnale, ma non in maniera drastica.
Allegato:
centeredDataPDF.jpg [ 23.87 KiB | Osservato 1296 volte ]
Guardando il valore stimato dal metodo si vede come questo pare comportarsi meglio del k-sigma. Si fa ingannare meno dagli outliers e non stima mai dati troppo distanti dal valore medio teorico di 125.
Pare essere un successo!
Allegato:
ExpectedValuePDF.jpg [ 20.63 KiB | Osservato 1296 volte ]
Nelle prossime puntate...Adesso sto provando ad applicare il metodo ai dati di Mauro che gentilmente mi ha fornito un centinaio di FITS. Innanzitutto voglio rispondere alle due domande sopra-citate.
L'applicazione dovrebbe essere abbastanza straight-forward. MA...!
Le immagini che ho purtroppo
non sono allineate. Questo fa sì che non abbia senso stackarle se prima non le allineo e...
non ci riesco! Sto provando a farlo con MATLAB, ma per ora non sono riuscito a cavare un ragno dal buco. Proverò ancora nei prossimi giorni, se però qualcuno ha qualche idea la ascolto molto volentieri.
A presto
Ciao
Luca