spectralpython / spectral

Python module for hyperspectral image processing
MIT License
573 stars 139 forks source link

Error at computing mnf #60

Closed aloboa closed 7 years ago

aloboa commented 7 years ago

I get an error at computing mnf. I've made test images in envi format by arranging 2 sets of 69 spectra as 2 images of 1 row x 69 columns x 50 bands. One image (testradVNIR) is for target radiance and the other (testdarkVNIR) for the dark noise. I have more bands but I have selected 50 to keep it smaller and avoid having more bands than spectra. Test images are available here https://dl.dropboxusercontent.com/u/3180464/test.zip

This is what I do:

from spectral import * data1 = envi.open('/media/alobo/LACIE500/Spectroradiometry/ADRIA/Tests/MNFdenoising/testradVNIR.hdr', '/media/alobo/LACIE500/Spectroradiometry/ADRIA/Tests/MNFdenoising/testradVNIR.envi').load() data2 = envi.open('/media/alobo/LACIE500/Spectroradiometry/ADRIA/Tests/MNFdenoising/testdarkVNIR.hdr', '/media/alobo/LACIE500/Spectroradiometry/ADRIA/Tests/MNFdenoising/testdarkVNIR.envi').load()

signal = calc_stats(data1) noise = calc_stats(data2) mnfr = mnf(signal, noise)

But get the following error

mnfr = mnf(signal, noise)

LinAlgError Traceback (most recent call last)

in () ----> 1 mnfr = mnf(signal, noise) /usr/local/lib/python2.7/dist-packages/spectral/algorithms/algorithms.pyc in mnf(signal, noise) 1644 from spectral.algorithms.algorithms import PrincipalComponents, GaussianStats 1645 C = noise.sqrt_inv_cov.dot(signal.cov).dot(noise.sqrt_inv_cov) -> 1646 (L, V) = np.linalg.eig(C) 1647 # numpy says eigenvalues may not be sorted so we'll sort them, if needed. 1648 if not np.alltrue(np.diff(L) <= 0): /usr/lib/python2.7/dist-packages/numpy/linalg/linalg.pyc in eig(a) 1016 _assertRank2(a) 1017 _assertSquareness(a) -> 1018 _assertFinite(a) 1019 a, t, result_t = _convertarray(a) # convert to double or cdouble type 1020 a = _to_native_byte_order(a) /usr/lib/python2.7/dist-packages/numpy/linalg/linalg.pyc in _assertFinite(*arrays) 163 for a in arrays: 164 if not (isfinite(a).all()): --> 165 raise LinAlgError("Array must not contain infs or NaNs") 166 167 def _assertNonEmpty(*arrays): LinAlgError: Array must not contain infs or NaNs I've made sure there are no infs or NaNs in my input data, it looks like a problem in the eigenanalysis. What can I do?
tboggs commented 7 years ago

It appears that either signal.cov or noise.cov may be singular. Even though you have more samples than bands, it may be that some of your spectra aren't linearly independent of the others. Off hand, I'd say try plotting their covariances (e.g. spy.imshow(signal.cov)) or checking to see that there isn't a NaN or inf in the covariance data. If those look fine try computing np.linalg.inv(signal.cov) and see what happens (and same for noise.cov). If computing the inverse fails you may have to either get more samples or regularize your data.

Thomas

On Wed, Nov 30, 2016 at 10:49 AM, aloboa notifications@github.com wrote:

I get an error at computing mnf. I've made test images in envi format by arranging 2 sets of 69 spectra as 2 images of 1 row x 69 columns x 50 bands. One image (testradVNIR) is for target radiance and the other (testdarkVNIR) for the dark noise. I have more bands but I have selected 50 to keep it smaller and avoid having more bands than spectra. Test images are available here https://dl.dropboxusercontent.com/u/3180464/test.zip

This is what I do:

from spectral import * data1 = envi.open('/media/alobo/LACIE500/Spectroradiometry/ ADRIA/Tests/MNFdenoising/testradVNIR.hdr', '/media/alobo/LACIE500/Spectroradiometry/ADRIA/Tests/ MNFdenoising/testradVNIR.envi').load() data2 = envi.open('/media/alobo/LACIE500/Spectroradiometry/ ADRIA/Tests/MNFdenoising/testdarkVNIR.hdr', '/media/alobo/LACIE500/Spectroradiometry/ADRIA/Tests/ MNFdenoising/testdarkVNIR.envi').load()

signal = calc_stats(data1) noise = calc_stats(data2) mnfr = mnf(signal, noise)

But get the following error mnfr = mnf(signal, noise)

LinAlgError Traceback (most recent call last) in () ----> 1 mnfr = mnf(signal, noise)

/usr/local/lib/python2.7/dist-packages/spectral/algorithms/algorithms.pyc in mnf(signal, noise) 1644 from spectral.algorithms.algorithms import PrincipalComponents, GaussianStats 1645 C = noise.sqrt_inv_cov.dot(signal.cov).dot(noise.sqrt_inv_cov) -> 1646 (L, V) = np.linalg.eig(C) 1647 # numpy says eigenvalues may not be sorted so we'll sort them, if needed. 1648 if not np.alltrue(np.diff(L) <= 0):

/usr/lib/python2.7/dist-packages/numpy/linalg/linalg.pyc in eig(a) 1016 _assertRank2(a) 1017 _assertSquareness(a) -> 1018 _assertFinite(a) 1019 a, t, result_t = _convertarray(a) # convert to double or cdouble type 1020 a = _to_native_byte_order(a)

/usr/lib/python2.7/dist-packages/numpy/linalg/linalg.pyc in _assertFinite(arrays) 163 for a in arrays: 164 if not (isfinite(a).all()): --> 165 raise LinAlgError("Array must not contain infs or NaNs") 166 167 def _assertNonEmpty(arrays):

LinAlgError: Array must not contain infs or NaNs

I've made sure there are no infs or NaNs in my input data, it looks like a problem in the eigenanalysis. What can I do?

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/spectralpython/spectral/issues/60, or mute the thread https://github.com/notifications/unsubscribe-auth/AEvuIZ7114LjFJfVfV9TAKc4yosaHjWhks5rDZr0gaJpZM4LAWF6 .

aloboa commented 7 years ago

Yes, the problem is on my side. A more clear error message in case of matrix singularity would be helpful though.

tboggs commented 7 years ago

I agree it would be useful to have a more informative error message but I would defer that to numpy where the eigenvalues are computed. I suppose an option is to catch LinAlgError but that is stated to be a generic exception class so I wouldn't want to mistakenly claim a singular matrix if in fact something else is occurring.

I'm going to close this for now but if someone comes up with a solution that doesn't introduce a significant performance hit, we can re-open it.