Load the Landsat OLI relative spectral responses
505874
def _load(self, scale=0.001):
        

        with open_workbook(self.path) as wb_:
            for sheet in wb_.sheets():
                if sheet.name in ['Plot of AllBands', ]:
                    continue
                ch_name = OLI_BAND_NAMES.get(sheet.name.strip())

                if ch_name != self.bandname:
                    continue

                wvl = sheet.col_values(0, 2)
                resp = sheet.col_values(1, 2)

                self.rsr = {'wavelength': np.array(wvl) / 1000.,
                            'response': np.array(resp)}
                break
610966
def _ir_calibrate(self, data, channel_name, cal_type):
        """"""
        if cal_type == 1:
            # spectral radiances
            return self._srads2bt(data, channel_name)
        elif cal_type == 2:
            # effective radiances
            return self._erads2bt(data, channel_name)
        else:
            raise NotImplementedError('Unknown calibration type')
944001
def load_ref_spectra():
    
    data_dir = "/Users/annaho/Data/AAOmega/ref_spectra"
    # Load the files & count the number of training objects
    ff = glob.glob("%s/*.txt" %data_dir)
    nstars = len(ff)
    print("We have %s training objects" %nstars)
    
    # Read the first file to get the wavelength array
    f = ff[0]
    data = Table.read(f, format="ascii.fast_no_header")
    wl = data['col1']
    npix = len(wl)
    print("We have %s pixels" %npix)

    tr_flux = np.zeros((nstars,npix))
    tr_ivar = np.zeros(tr_flux.shape)

    for i,f in enumerate(ff):
        data = Table.read(f, format="ascii.fast_no_header")
        flux = data['col2']
        tr_flux[i,:] = flux
        sigma = data['col3']
        tr_ivar[i,:] = 1.0 / sigma**2

    return np.array(ff), wl, tr_flux, tr_ivar
964188
def ndvi(self):
        
        red, nir = self.reflectance(3), self.reflectance(4)
        ndvi = self._divide_zero((nir - red), (nir + red), nan)

        return ndvi
122622
private void thresholdBlur(Kernel kernel, int[] inPixels, int[] outPixels, int width, int height, boolean alpha) {
		int index = 0;
		float[] matrix = kernel.getKernelData( null );
		int cols = kernel.getWidth();
		int cols2 = cols/2;

		for (int y = 0; y < height; y++) {
			int ioffset = y*width;
            int outIndex = y;
			for (int x = 0; x < width; x++) {
				float r = 0, g = 0, b = 0, a = 0;
				int moffset = cols2;

                int rgb1 = inPixels[ioffset+x];
                int a1 = (rgb1 >> 24) & 0xff;
                int r1 = (rgb1 >> 16) & 0xff;
                int g1 = (rgb1 >> 8) & 0xff;
                int b1 = rgb1 & 0xff;
				float af = 0, rf = 0, gf = 0, bf = 0;
                for (int col = -cols2; col <= cols2; col++) {
					float f = matrix[moffset+col];

					if (f != 0) {
						int ix = x+col;
						if (!(0 <= ix && ix < width))
							ix = x;
						int rgb2 = inPixels[ioffset+ix];
                        int a2 = (rgb2 >> 24) & 0xff;
                        int r2 = (rgb2 >> 16) & 0xff;
                        int g2 = (rgb2 >> 8) & 0xff;
                        int b2 = rgb2 & 0xff;

						int d;
                        d = a1-a2;
                        if ( d >= -threshold && d <= threshold ) {
                            a += f * a2;
                            af += f;
                        }
                        d = r1-r2;
                        if ( d >= -threshold && d <= threshold ) {
                            r += f * r2;
                            rf += f;
                        }
                        d = g1-g2;
                        if ( d >= -threshold && d <= threshold ) {
                            g += f * g2;
                            gf += f;
                        }
                        d = b1-b2;
                        if ( d >= -threshold && d <= threshold ) {
                            b += f * b2;
                            bf += f;
                        }
					}
				}
                a = af == 0 ? a1 : a/af;
                r = rf == 0 ? r1 : r/rf;
                g = gf == 0 ? g1 : g/gf;
                b = bf == 0 ? b1 : b/bf;
				int ia = alpha ? PixelUtils.clamp((int)(a+0.5)) : 0xff;
				int ir = PixelUtils.clamp((int)(r+0.5));
				int ig = PixelUtils.clamp((int)(g+0.5));
				int ib = PixelUtils.clamp((int)(b+0.5));
				outPixels[outIndex] = (ia << 24) | (ir << 16) | (ig << 8) | ib;
                outIndex += height;
			}
		}
	}
358036
def _load(self, scale=1.0):
        

        LOG.debug("File: %s", str(self.requested_band_filename))
        ncf = Dataset(self.requested_band_filename, 'r')

        wvl = ncf.variables['wavelength'][:] * scale
        resp = ncf.variables['response'][:]

        self.rsr = {'wavelength': wvl, 'response': resp}
1274729
def _load_calib(self):
        """"""
        # We'll build the calibration parameters as a dictionary, then
        # convert it to a namedtuple to prevent it from being modified later
        data = {}

        # Load the rigid transformation from IMU to velodyne
        data['T_velo_imu'] = self._load_calib_rigid('calib_imu_to_velo.txt')

        # Load the camera intrinsics and extrinsics
        data.update(self._load_calib_cam_to_cam(
            'calib_velo_to_cam.txt', 'calib_cam_to_cam.txt'))

        # Pre-compute the IMU to rectified camera coordinate transforms
        data['T_cam0_imu'] = data['T_cam0_velo'].dot(data['T_velo_imu'])
        data['T_cam1_imu'] = data['T_cam1_velo'].dot(data['T_velo_imu'])
        data['T_cam2_imu'] = data['T_cam2_velo'].dot(data['T_velo_imu'])
        data['T_cam3_imu'] = data['T_cam3_velo'].dot(data['T_velo_imu'])

        self.calib = namedtuple('CalibData', data.keys())(*data.values())
1064814
def _load(self, scale=1.0):
        """"""
        data = np.genfromtxt(self.requested_band_filename,
                             unpack=True,
                             names=['wavenumber',
                                    'response'],
                             skip_header=4)

        # Data are wavenumbers in cm-1:
        wavelength = 1. / data['wavenumber'] * 10000.
        response = data['response']

        # The real MetImage has 24 detectors. However, for now we store the
        # single rsr as 'detector-1', indicating that there will be multiple
        # detectors in the future:
        detectors = {}
        detectors['det-1'] = {'wavelength': wavelength, 'response': response}
        self.rsr = detectors
935828
def get_reflectance(self, sun_zenith, sat_zenith, azidiff, bandname, redband=None):
        """"""
        # Get wavelength in nm for band:
        if isinstance(bandname, float):
            LOG.warning('A wavelength is provided instead of band name - ' +
                        'disregard the relative spectral responses and assume ' +
                        'it is the effective wavelength: %f (micro meter)', bandname)
            wvl = bandname * 1000.0
        else:
            wvl = self.get_effective_wavelength(bandname)
            wvl = wvl * 1000.0

        rayl, wvl_coord, azid_coord, satz_sec_coord, sunz_sec_coord = self.get_reflectance_lut()

        # force dask arrays
        compute = False
        if HAVE_DASK and not isinstance(sun_zenith, Array):
            compute = True
            sun_zenith = from_array(sun_zenith, chunks=sun_zenith.shape)
            sat_zenith = from_array(sat_zenith, chunks=sat_zenith.shape)
            azidiff = from_array(azidiff, chunks=azidiff.shape)
            if redband is not None:
                redband = from_array(redband, chunks=redband.shape)

        clip_angle = rad2deg(arccos(1. / sunz_sec_coord.max()))
        sun_zenith = clip(sun_zenith, 0, clip_angle)
        sunzsec = 1. / cos(deg2rad(sun_zenith))
        clip_angle = rad2deg(arccos(1. / satz_sec_coord.max()))
        sat_zenith = clip(sat_zenith, 0, clip_angle)
        satzsec = 1. / cos(deg2rad(sat_zenith))
        shape = sun_zenith.shape

        if not(wvl_coord.min() < wvl < wvl_coord.max()):
            LOG.warning(
                "Effective wavelength for band %s outside 400-800 nm range!",
                str(bandname))
            LOG.info(
                "Set the rayleigh/aerosol reflectance contribution to zero!")
            if HAVE_DASK:
                chunks = sun_zenith.chunks if redband is None else redband.chunks
                res = zeros(shape, chunks=chunks)
                return res.compute() if compute else res
            else:
                return zeros(shape)

        idx = np.searchsorted(wvl_coord, wvl)
        wvl1 = wvl_coord[idx - 1]
        wvl2 = wvl_coord[idx]

        fac = (wvl2 - wvl) / (wvl2 - wvl1)
        raylwvl = fac * rayl[idx - 1, :, :, :] + (1 - fac) * rayl[idx, :, :, :]
        tic = time.time()

        smin = [sunz_sec_coord[0], azid_coord[0], satz_sec_coord[0]]
        smax = [sunz_sec_coord[-1], azid_coord[-1], satz_sec_coord[-1]]
        orders = [
            len(sunz_sec_coord), len(azid_coord), len(satz_sec_coord)]
        f_3d_grid = atleast_2d(raylwvl.ravel())

        if HAVE_DASK and isinstance(smin[0], Array):
            # compute all of these at the same time before passing to the interpolator
            # otherwise they are computed separately
            smin, smax, orders, f_3d_grid = da.compute(smin, smax, orders, f_3d_grid)
        minterp = MultilinearInterpolator(smin, smax, orders)
        minterp.set_values(f_3d_grid)

        if HAVE_DASK:
            ipn = map_blocks(self._do_interp, minterp, sunzsec, azidiff,
                             satzsec, dtype=raylwvl.dtype, chunks=azidiff.chunks)
        else:
            ipn = self._do_interp(minterp, sunzsec, azidiff, satzsec)

        LOG.debug("Time - Interpolation: {0:f}".format(time.time() - tic))

        ipn *= 100
        res = ipn
        if redband is not None:
            res = where(redband < 20., res,
                        (1 - (redband - 20) / 80) * res)

        res = clip(res, 0, 100)
        if compute:
            res = res.compute()
        return res
332897
def stand_alone_imagery(self):
        
        return json.loads(self._query(IMAGE, FORECAST, SURFACE_PRESSURE, CAPABILITIES, "").decode(errors="replace"))
1282382
function fftImage(image) {
    const real = createMatrix();
    const imag = createMatrix();
    let idx = 0;
    const nchan = image.channels;
    const alphachan = 1 - (nchan % 2);
    const colorchan = nchan - alphachan;
    for (let y = 0; y < image.height; ++y) {
        for (let x = 0; x < image.width; ++x) {
            let v = 0;
            for (let c = 0; c < colorchan; ++c) {
                v += 255 - image.data[idx++];
            }
            for (let c = 0; c < alphachan; ++c) {
                v += image.data[idx++];
            }
            real.set(y, x, v);
        }
    }
    fft(1, real, imag);
    return {
        real: real,
        imag: imag,
        width: image.width,
        height: image.height,
    };
}
1250196
def _vis_calibrate(self, data, key):
        """"""
        # radiance to reflectance taken as in mipp/xrit/MSG.py
        # again FCI User Guide is not clear on how to do this

        sirr = self.nc[
            '/data/{}/measured/channel_effective_solar_irradiance'
            .format(key.name)][...]

        # reflectance = radiance / sirr * 100
        data.data[:] /= sirr
        data.data[:] *= 100
645641
def _load(self, scale=1.0):
        

        LOG.debug("File: %s", str(self.requested_band_filename))

        data = np.genfromtxt(self.requested_band_filename,
                             unpack=True,
                             names=['wavelength',
                                    'wavenumber',
                                    'response'],
                             skip_header=2)

        wvl = data['wavelength'] * scale
        resp = data['response']

        self.rsr = {'wavelength': wvl, 'response': resp}
107637
def make_tb2rad_lut(self, filepath, normalized=True):
        """"""
        tb_ = np.arange(TB_MIN, TB_MAX, self.tb_resolution)
        retv = self.tb2radiance(tb_, normalized=normalized)
        rad = retv['radiance']
        np.savez(filepath, tb=tb_, radiance=rad)
255315
def _get_binned_arrays(self, wavelengths, flux_unit, area=None,
                           vegaspec=None):
        """"""
        x = self._validate_binned_wavelengths(wavelengths)
        y = self.sample_binned(wavelengths=x, flux_unit=flux_unit, area=area,
                               vegaspec=vegaspec)

        if isinstance(wavelengths, u.Quantity):
            w = x.to(wavelengths.unit, u.spectral())
        else:
            w = x

        return w, y
15676
def curvature(self, name, **kwargs):
        

        name = self.roi.get_source_by_name(name).name

        saved_state = LikelihoodState(self.like)
        source = self.components[0].like.logLike.getSource(str(name))

        old_spectrum = source.spectrum()
        old_pars = copy.deepcopy(self.roi[name].spectral_pars)
        old_type = self.roi[name]['SpectrumType']

        if old_type != 'PowerLaw':

            dnde = self.like[name].spectrum()(pyLike.dArg(1000.))
            value, scale = utils.scale_parameter(dnde)
            pars0 = {
                'Prefactor':
                    {'value': value, 'scale': scale,
                        'min': 1E-5, 'max': 1000., 'free': True},
                'Index':
                    {'value': 2.0, 'scale': -1.0, 'min': 0.0,
                        'max': 5.0, 'free': False},
                'Scale':
                    {'value': 1E3, 'scale': 1.0, 'min': 1.,
                        'max': 1E6, 'free': False},
            }
            self.set_source_spectrum(str(name), 'PowerLaw',
                                     spectrum_pars=pars0,
                                     update_source=False)

        self.free_source(name, loglevel=logging.DEBUG)
        fit_pl = self._fit(loglevel=logging.DEBUG)

        prefactor = self._get_param(name, 'Prefactor')
        index = self._get_param(name, 'Index')
        scale = self._get_param(name, 'Scale')

        pars1 = {
            'norm': copy.deepcopy(prefactor),
            'alpha': copy.deepcopy(index),
            'Eb': copy.deepcopy(scale),
        }
        pars1['alpha']['scale'] *= -1
        pars1['alpha']['min'] = -5.0
        pars1['alpha']['max'] = 5.0

        self.set_source_spectrum(str(name), 'LogParabola',
                                 spectrum_pars=pars1,
                                 update_source=False)

        self.free_source(name, loglevel=logging.DEBUG)
        fit_lp = self._fit(loglevel=logging.DEBUG)

        self.set_source_spectrum(str(name), old_type,
                                 spectrum_pars=old_pars,
                                 update_source=False)

        pars2 = {
            'Prefactor': copy.deepcopy(prefactor),
            'Index1': copy.deepcopy(index),
            'Cutoff': {'value': 1000.0, 'scale': 1E3,
                       'min': 10.0, 'max': 1E4, 'free': True},
            'Index2': {'value': 1.0, 'scale': 1.0,
                       'min': 1.0, 'max': 1.0, 'free': False},
            'Scale': copy.deepcopy(scale)
        }

        self.set_source_spectrum(str(name), 'PLSuperExpCutoff',
                                 spectrum_pars=pars2,
                                 update_source=False)

        self.free_source(name, loglevel=logging.DEBUG)
        fit_ple = self._fit(loglevel=logging.DEBUG)

        # Revert to initial spectral model
        self.set_source_spectrum(str(name), old_type,
                                 spectrum_pars=old_pars,
                                 update_source=False)
        saved_state.restore()

        lp_ts_curv = 2.0 * (fit_lp['loglike'] - fit_pl['loglike'])
        ple_ts_curv = 2.0 * (fit_ple['loglike'] - fit_pl['loglike'])
        o = MutableNamedTuple(ts_curv=lp_ts_curv,
                              lp_ts_curv=lp_ts_curv,
                              ple_ts_curv=ple_ts_curv,
                              loglike_pl=fit_pl['loglike'],
                              loglike_lp=fit_lp['loglike'],
                              loglike_ple=fit_ple['loglike'])

        self.logger.info('LogLike_PL: %12.3f LogLike_LP: %12.3f LogLike_PLE: %12.3f',
                         o.loglike_pl, o.loglike_lp, o.loglike_ple)
        self.logger.info('TS_curv:        %.3f (LP)', o.lp_ts_curv)
        self.logger.info('TS_curv:        %.3f (PLE)', o.ple_ts_curv)
        return o
893125
private void computeOctaveScales() {
		octaveImages[0] = tempImage0;
		for (int i = 1; i < numScales+3; i++) {
			octaveImages[i].reshape(tempImage0.width, tempImage0.height);
			applyGaussian(octaveImages[i - 1], octaveImages[i], kernelSigmaToK[i-1]);
		}

		for (int i = 1; i < numScales+3; i++) {
			differenceOfGaussian[i-1].reshape(tempImage0.width, tempImage0.height);
			PixelMath.subtract(octaveImages[i],octaveImages[i - 1],differenceOfGaussian[i-1]);
		}
	}