Sto usando numpy e pyfits per manipolare gli spettri e ho bisogno di alta precisione (qualcosa come 8-10 posizioni decimali su un valore che potrebbe arrivare fino a 10^12). Per questo il tipo di dati "decimali" sarebbe perfetto (float64 non è abbastanza buono), ma numpy.interp unfortunalely non piace:Numpy alta precisione
File ".../modules/manip_fits.py", line 47, in get_shift
pix_shift = np.interp(x, xp, fp)-fp
File "/usr/lib/python2.7/dist-packages/numpy/lib/function_base.py", line 1053, in interp
return compiled_interp(x, xp, fp, left, right)
TypeError: array cannot be safely cast to required type
Una versione semplificata del codice che uso:
fp = np.array(range(new_wave.shape[-1]),dtype=Decimal)
pix_shift = np.empty_like(wave,dtype=Decimal)
x = wave
xp = new_wave
pix_shift = np.interp(x, xp, fp)-fp
dove 'wave' e 'new_wave' sono una matrice numpy a una dimensione che rappresenta uno spettro 1D. Questo codice è necessario per spostare i miei spettri lungo l'asse x (che è la lunghezza d'onda)
Il mio problema più grande è che più in basso nel codice divido i miei spettri con uno spettro modello costruito dalla somma di tutti i miei spettri per analizzare le differenze e poiché non ho abbastanza cifre decimali ho degli errori di arrotondamento. Qualche idea?
Grazie!
UPDATE:
test Esempio:
import numpy as np
from decimal import *
getcontext().prec = 12
wave = np.array([Decimal(xx*np.pi) for xx in range(0,10)],dtype=np.dtype(Decimal))
new_wave = np.array([Decimal(xx*np.pi+0.5) for xx in range(0,10)],dtype=np.dtype(Decimal))
fp = np.array(range(new_wave.shape[-1]),dtype=Decimal)
pix_shift = np.empty_like(wave,dtype=Decimal)
x = wave
xp = new_wave
pix_shift = np.interp(x, xp, fp)-fp
L'errore è:
Traceback (most recent call last):
File "untitled.py", line 16, in <module>
pix_shift = np.interp(x, xp, fp)-fp
File "/usr/lib/python2.7/dist-packages/numpy/lib/function_base.py", line 1053, in interp
return compiled_interp(x, xp, fp, left, right)
TypeError: array cannot be safely cast to required type
questo è il più vicino posso fornire senza utilizzare il vero spettri in formato FITS.
UPDATE 2: Alcuni valori tipici della mia spettri, stampata utilizzando decimale:
18786960689.118938446044921875
18473926205.282184600830078125
18325454516.792461395263671875
18400241010.149127960205078125
2577901751996.03857421875
2571812230557.63330078125
2567431795280.80712890625
il problema che sto ottenendo è quando faccio le operazioni tra di loro, vengo arrotondamenti errori. Ad esempio, creo un modello per tutti gli spettri sommandoli tutti. Quindi uso questo modello per normalizzare ogni spettro. Un esempio:
Spectra = np.array([Spectrum1, Spectrum2, ...])
Template = np.nansum(Spectra, axis= 0)
NormSpectra = Spectra/Template
Questo mi deve restituire solo il rumore del spettri (supponendo che il modello è una buona rappresentazione della stella). Ho provato normalizzare ogni spettro al suo flusso totale
(Spectrum1 = Spectrum1/np.nansum(Spectrum1), ...)
così come template, ma sarebbe ottenere ancora peggiori errori di arrotondamento in su.
L'utilizzo di Decimal funzionerebbe correttamente per me, ma ho bisogno di "spostare" i miei spettri in modo che tutte le caratteristiche/linee spettrali siano allineate.
Spero che questo abbia senso?
Hai provato a usare 'numpy.longdouble' come dtype? http://mail.scipy.org/pipermail/scipy-dev/2008-March/008562.html –
@Zhenya: ho appena provato, e sto facendo lo stilin "TypeError: array non può essere lanciato in modo sicuro nel tipo richiesto" :( – jorgehumberto
puoi fornire un piccolo esempio autosufficiente che dimostri il problema? –