how to indentify infinite values in array

[bobc]

Hi Joleen -

I'm having some trouble replicating your problem from your 02/10 post below:

This range of this ndvi should be -1 to 1. However, you see that the range is off by an order of magnitude. However, if the image is displayed and the data values are probed, it looks like the data values in the display are correct.

Here's the script I'm running:

Code: Select all

parmsList=dict(
    server='dcarchivex.ssec.wisc.edu',
    dataset='GREX',
    descriptor='FD',
    band=1,
    day='2016277',
    time=('16:04', '16:09'),
    position='ALL'
    )

dateTimes=listADDEImageTimes(**parmsList)

if (len(dateTimes) >= 1):
    #panel=buildWindow(height=900,width=900)[0]

    for dateTime in dateTimes[::2]:
       index=1
       parms=dict(
           server=parmsList['server'],
           dataset=parmsList['dataset'],
           descriptor='FD',
           day=dateTime['day'],
           time=dateTime['time'],
           unit='ALB',
           location=(0.1,-98.4),
           coordinateSystem=LATLON,
           place=ULEFT,
           mag=(1,1)
       )
       rdata=loadADDEImage(band=2, size=(8,8), **parms)
       gdata=loadADDEImage(band=3, size=(4,4), **parms)
       bdata=loadADDEImage(band=4, size=(2,2), **parms)
       
       rfield=getRangeType(rdata)
       bfield=getRangeType(bdata)
       gfield=getRangeType(gdata)
       
       #rUnit=(rfield.defaultUnits)[0]
       #bUnit=(bfield.defaultUnits)[0]
       #nIRUnit=(nIRfield.defaultUnits)[0]
       rUnit=whatType(rdata).toString().rstrip().splitlines()[-1].split('Unit: ',1)[1]
       gUnit=whatType(gdata).toString().rstrip().splitlines()[-1].split('Unit: ',1)[1]
       bUnit=whatType(bdata).toString().rstrip().splitlines()[-1].split('Unit: ',1)[1]
       
       rdata=newUnit(rdata,rfield.toString(),None)
       bdata=newUnit(bdata,bfield.toString(),None)
       gdata=newUnit(gdata,gfield.toString(),None)

       ndviDenom = sub(gdata,-rdata)
       ndviNumer = sub(gdata,rdata)
       ndviData = ndviNumer/ndviDenom

ndviLayer = activeDisplay().createLayer('Image Display',ndviData)
gridLayer = activeDisplay().createLayer('Grid Table', ndviData)

Notes:

• I changed nIRfield/nIRdata to 'gfield/gdata'.
• I hadn't used the bluefield/bluedata earlier on so I added them to this script, though I'm not actually using them in the NDVI calculation.

Once I'm done the script, the enhancement range in the Main Display goes from -1 to 0.13, which matches the Grid Table output. I also ran 'describe(ndviData)' which had these same min/max values. I'm wondering if you are doing something extra in your script than I am, which looks to be possible since you are working with the bfield/bdata objects in the bit of code from your last post.

Thanks -
Bob

[joleenf]

Hi Bob,

It will take me a bit to get back to this. We are opting to go with another routine to achieve similar results. However, when I get a chance, I will try to recreate the problem on my machine. I never saw this problem with DOE-4 test data, which is what you are running in the script below. In fact, the day/time you were running, provided very nice results for me. I saw this problem with checkout data from GOES-16 and in some AHI sample files from Ray.

I will get back to you on a day/time combo fro G16.

Joleen

[joleenf]

Hi Bob,

I apologize, I could not get back far enough in my history to find the GOES16 data that was giving me problems. However, I think you have solved this issue for me and there are a few things to note in summary:

1.) If I am controlling a function, such as NDVI, I should mask the denominator before dividing (makes sense, should have done that first).

2.) If I don't have full control over the function, I can still use the formula which is creating infinite values, and if I want to mask the infinite values to missing, I can divide the result by the abs_data(result) to create a mask of missing and 1.

ndvi=ndvi*(ndvi)/abs_data(ndvi)

**abs_data should be used so that it does no conflict with Jython built-in.

3.) Be careful when using rescale. As noted in the documentation:

f - the FieldImpl or FlatField
inlo - the input low-range value
inhi - the input high-range value
outlo - the output low-range value
outhi - the output high range value Values of the original field will be linearly scaled from "inlo:inhi" to "outlo:outhi" Values < inlo will be set to outlo; values > inhi set to outhi If input FieldImpl is a sequence, then all items in sequence are done

This means that -infinity is set to outlo, and infinity is set to outhi.

Also in rescale, outlo must be less than outhi. Otherwise, all values are set to zero...

Examples:

Code: Select all

import java.lang.Float

infinity=Float.POSITIVE_INFINITY
neg_infinity=Float.NEGATIVE_INFINITY

a=field([infinity, neg_infinity, 1, 0, 2])


strA='\n'.join('Original Field: {}'.format(k) for i,k in enumerate(a.getValues()))

print strA

b=rescale(a,-1,1,0,10)

strb='\n'.join('Rescale Field -1,1 to 0,10: {}'.format(k) for i,k in enumerate(b.getValues()))

print strb

c=rescale(a,-1,1,10,0)

strc='\n'.join('Rescale Field -1,1 to 10,0: {}'.format(k) for i,k in enumerate(c.getValues()))

print 'THIS DOES NOT WORK!'
print strc

d=rescale(a,1,-1,0,10)

strd='\n'.join('Rescale Field 1,-1 to 0,10: {}'.format(k) for i,k in enumerate(d.getValues()))

print 'Reverse inlo, inhi to flip low and hi in output field'
print strd

e=rescale(a,1,-1,10,0)

stre='\n'.join('Rescale Field 1,-1 to 10,0: {}'.format(k) for i,k in enumerate(e.getValues()))

print stre

Thanks,
Joleen