2013-05-15 3 views
9

나는 numpy와 pyfits를 사용하여 스펙트럼을 조작하고 고정밀도 (10^12만큼 높을 수있는 값에 8-10 소수 자릿수가 필요함)가 필요합니다.Numpy 고정밀

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 

코드의 단순화 된 버전은 내가 사용 :

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 
그 데이터 유형 "소수"에 대한 (float64가 충분하지 않습니다)하지만, unfortunalely numpy.interp이 그것을 좋아하지 않는 완벽 할 것

여기서 'wave'및 'new_wave'는 1D 스펙트럼을 나타내는 일차원 numpy 배열입니다. 이 코드는 x 축을 따라 내 스펙트럼을 이동하는 데 필요합니다 (wavelenght)

내 가장 큰 문제는 코드를 더 내려가는 것입니다. 모든 스펙트럼의 합계로 구성된 템플릿 스펙트럼으로 내 스펙트럼을 나눕니다. 차이점을 분석하고 소수점이 충분하지 않아 반올림 오류가 발생합니다. 어떤 아이디어?

감사합니다.

UPDATE :

시험 예 :

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 

오류는 다음과 같습니다

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 

이것은 내가에 실제 스펙트럼을 사용하지 않고 제공 할 수있는 가장 가까운 형식을 맞는다.

업데이트 2 : 내 스펙트럼의 몇 가지 전형적인 값, 진수 사용하여 인쇄 : 나는 그들 사이의 작업, 내가 오류를 반올림 취득을 할 때

18786960689.118938446044921875 
    18473926205.282184600830078125 
    18325454516.792461395263671875 
    18400241010.149127960205078125 
2577901751996.03857421875 
2571812230557.63330078125 
2567431795280.80712890625 

나는 점점 오전 문제입니다. 예를 들어 모든 스펙트럼을 합산하여 모든 스펙트럼에 대한 템플릿을 만듭니다. 그런 다음이 템플릿을 사용하여 모든 스펙트럼을 정규화합니다. 예 :

Spectra = np.array([Spectrum1, Spectrum2, ...]) 
Template = np.nansum(Spectra, axis= 0) 

NormSpectra = Spectra/Template 

이 (템플릿이 별의 좋은 표현이라고 가정) 스펙트럼에 나에게 단지 소음을 반환해야합니다. 각 스펙트럼을 총 플럭스로 표준화하려고 시도했지만 템플릿은 물론 템플릿을 더 잘 반올림하여 오류를 반올림합니다.

Decimal을 사용하면 나에게 잘 작동하지만 모든 스펙트럼 피처/선이 정렬되도록 내 스펙트럼을 "이동"해야합니다.

희망이 있으십니까?

+0

'numpy.longdouble'을 dtype으로 사용해 보셨습니까? http://mail.scipy.org/pipermail/scipy-dev/2008-March/008562.html –

+0

@Zhenya : 방금 시도해 보았습니다. "TypeError : 배열을 필수 형식으로 안전하게 캐스팅 할 수 없습니다." – jorgehumberto

+0

문제를 설명하는 작은 자체 포함 예제를 제공 할 수 있습니까? –

답변

5

np.float64을 어떻게 확인할 수 있습니까? 전형적인 유스 케이스에서는 이중에서 15 개의 유효 숫자를 기대할 수 있습니다.

이 정도면 충분하지 않다면 (별칭 np.longdouble)으로 시도해 볼 수 있습니다.

그러나 귀하의 문제는 그보다 더 심해 보입니다. 부적절한 문제 (작은 숫자로 큰 숫자를 나누는 것이 일반적 임) 인 것 같습니다. 그리고 이것은 당신이 원하는 어떤 것이 아닙니다.정밀도를 높이면 문제가 어느 정도 해결되지만 float256/float512/등이 필요한 데이터가 발생합니다. 병리학적인 반올림 오류를 방지합니다.

모든 문제 (XY Problem)를 해결할 수있는 또 다른 방법을 찾을 수 있도록 해결책이 아니라 문제를 설명하는 것이 좋습니다.

관련 문제