Return the complex amplitudes coupled into a set of fibers
195385
def FibreCouple(pupils,modeDiameter):
    
    gridSize=pupils.shape[-1]
    pupilsVector=np.reshape(pupils,(-1,gridSize**2))
    mode=np.reshape(FibreMode(gridSize,modeDiameter),(gridSize**2,))
    return np.inner(pupilsVector,mode)
978048
public static void multiplyRealComplex( GrayF realA ,
											ImageInterleaved complexB , ImageInterleaved complexC ) {
		if( complexB instanceof InterleavedF32 ) {
			DiscreteFourierTransformOps.multiplyRealComplex((GrayF32) realA, (InterleavedF32) complexB, (InterleavedF32) complexC);
		} else if( complexB instanceof InterleavedF64 ) {
			DiscreteFourierTransformOps.multiplyRealComplex((GrayF64) realA, (InterleavedF64) complexB, (InterleavedF64) complexC);
		} else {
			throw new IllegalArgumentException("Unknown image type");
		}
	}
916154
def choi(self) -> bk.BKTensor:
        """"""
        # Put superop axes in [ok, ib, ob, ik] and reshape to matrix
        N = self.qubit_nb
        return bk.reshape(self.sharp.tensor, [2**(N*2)] * 2)
1269882
def _coupling_matrix(self, lmax, nwin=None, weights=None):
        """"""
        if nwin is None:
            nwin = self.nwin

        if weights is None:
            weights = self.weights

        tapers_power = _np.zeros((self.lwin+1, nwin))
        for i in range(nwin):
            tapers_power[:, i] = _spectrum(self.to_array(i),
                                           normalization='4pi',
                                           convention='power', unit='per_l')

        if weights is None:
            return _shtools.SHMTCouplingMatrix(lmax, tapers_power, k=nwin)
        else:
            return _shtools.SHMTCouplingMatrix(lmax, tapers_power, k=nwin,
                                               taper_wt=self.weights)
295837
public static void multiplyComplex( InterleavedF32 complexA , InterleavedF32 complexB , InterleavedF32 complexC ) {

		InputSanityCheck.checkSameShape(complexA, complexB,complexC);

		for( int y = 0; y < complexA.height; y++ ) {

			int indexA = complexA.startIndex + y*complexA.stride;
			int indexB = complexB.startIndex + y*complexB.stride;
			int indexC = complexC.startIndex + y*complexC.stride;

			for( int x = 0; x < complexA.width; x++, indexA += 2 , indexB += 2  ,indexC += 2 ) {

				float realA = complexA.data[indexA];
				float imgA = complexA.data[indexA+1];
				float realB = complexB.data[indexB];
				float imgB = complexB.data[indexB+1];

				complexC.data[indexC] = realA*realB - imgA*imgB;
				complexC.data[indexC+1] = realA*imgB + imgA*realB;
			}
		}
	}
32165
def xstep(self):
        """

        self.YU[:] = self.Y - self.U
        YUf = sl.rfftn(self.YU, None, self.cri.axisN)
        YUf0 = self.block_sep0(YUf)
        YUf1 = self.block_sep1(YUf)

        b = self.rho * np.sum(np.conj(self.GDf) * YUf1, axis=-1)
        if self.cri.Cd > 1:
            b = np.sum(b, axis=self.cri.axisC, keepdims=True)
        b += self.DSf + self.rho*YUf0

        # Concatenate multiple GDf components on axisC. For
        # single-channel signals, and multi-channel signals with a
        # single-channel dictionary, we end up with sl.solvemdbi_ism
        # solving a linear system of rank dimN+1 (corresponding to the
        # dictionary and a gradient operator per spatial dimension) plus
        # an identity. For multi-channel signals with a multi-channel
        # dictionary, we end up with sl.solvemdbi_ism solving a linear
        # system of rank C.d (dimN+1) (corresponding to the dictionary
        # and a gradient operator per spatial dimension for each
        # channel) plus an identity.

        # The structure of the linear system to be solved depends on the
        # number of channels in the signal and dictionary. Both branches are
        # the same in the single-channel signal case (the choice of handling
        # it via the 'else' branch is somewhat arbitrary).
        if self.cri.C > 1 and self.cri.Cd == 1:
            # Concatenate multiple GDf components on the final axis
            # of GDf (that indexes the number of gradient operators). For
            # multi-channel signals with a single-channel dictionary,
            # sl.solvemdbi_ism has to solve a linear system of rank dimN+1
            # (corresponding to the dictionary and a gradient operator per
            # spatial dimension)
            DfGDf = np.concatenate(
                [self.Df[..., np.newaxis],] +
                [np.sqrt(self.rho)*self.GDf[..., k, np.newaxis] for k
                 in range(self.GDf.shape[-1])], axis=-1)
            self.Xf[:] = sl.solvemdbi_ism(DfGDf, self.rho, b[..., np.newaxis],
                                          self.cri.axisM, -1)[..., 0]
        else:
            # Concatenate multiple GDf components on axisC. For multi-channel
            # signals with a multi-channel dictionary, sl.solvemdbi_ism has
            # to solve a linear system of rank C.d (dimN+1) (corresponding to
            # the dictionary and a gradient operator per spatial dimension
            # for each channel) plus an identity.
            DfGDf = np.concatenate(
                [self.Df,] + [np.sqrt(self.rho)*self.GDf[..., k] for k
                              in range(self.GDf.shape[-1])],
                axis=self.cri.axisC)
            self.Xf[:] = sl.solvemdbi_ism(DfGDf, self.rho, b, self.cri.axisM,
                                          self.cri.axisC)

        self.X = sl.irfftn(self.Xf, self.cri.Nv, self.cri.axisN)

        if self.opt['LinSolveCheck']:
            if self.cri.C > 1 and self.cri.Cd == 1:
                Dop = lambda x: sl.inner(DfGDf, x[..., np.newaxis],
                                         axis=self.cri.axisM)
                DHop = lambda x: sl.inner(np.conj(DfGDf), x, axis=-1)
                ax = DHop(Dop(self.Xf))[..., 0] + self.rho*self.Xf
            else:
                Dop = lambda x: sl.inner(DfGDf, x, axis=self.cri.axisM)
                DHop = lambda x: sl.inner(np.conj(DfGDf), x,
                                          axis=self.cri.axisC)
                ax = DHop(Dop(self.Xf)) + self.rho*self.Xf
            self.xrrs = sl.rrs(ax, b)
        else:
            self.xrrs = None
854125
def _extract_intensities(image, mask = slice(None)):
    
    return numpy.array(image, copy=True)[mask].ravel()
1318912
def _series_sis(self):
        """"""
        pages = self.pages._getlist(validate=False)
        page = pages[0]
        lenpages = len(pages)
        md = self.sis_metadata
        if 'shape' in md and 'axes' in md:
            shape = md['shape'] + page.shape
            axes = md['axes'] + page.axes
        elif lenpages == 1:
            shape = page.shape
            axes = page.axes
        else:
            shape = (lenpages,) + page.shape
            axes = 'I' + page.axes
        self.is_uniform = True
        return [TiffPageSeries(pages, shape, page.dtype, axes, kind='SIS')]
1164261
def mixing_phases(U):
    """"""
    f = {}
    # angles
    f['t13'] = asin(abs(U[0,2]))
    if U[0,0] == 0:
        f['t12'] = pi/2
    else:
        f['t12'] = atan(abs(U[0,1])/abs(U[0,0]))
    if U[2,2] == 0:
        f['t23'] = pi/2
    else:
        f['t23'] = atan(abs(U[1,2])/abs(U[2,2]))
    s12 = sin(f['t12'])
    c12 = cos(f['t12'])
    s13 = sin(f['t13'])
    c13 = cos(f['t13'])
    s23 = sin(f['t23'])
    c23 = cos(f['t23'])
    # standard phase
    if (s12*s23) == 0 or (c12*c13**2*c23*s13) == 0:
        f['delta'] = 0
    else:
        f['delta'] = -phase((U[0,0].conj()*U[0,2]*U[2,0]*U[2,2].conj()/(c12*c13**2*c23*s13) + c12*c23*s13)/(s12*s23))
    # Majorana phases
    f['delta1']  = phase(exp(1j*f['delta']) * U[0, 2])
    f['delta2']  = phase(U[1, 2])
    f['delta3']  = phase(U[2, 2])
    f['phi1'] = 2*phase(exp(1j*f['delta1']) * U[0, 0].conj())
    f['phi2'] = 2*phase(exp(1j*f['delta1']) * U[0, 1].conj())
    return f
61712
def blend(c1, c2):
    """"""
    return [c1[i] * (0xFF - c2[3]) + c2[i] * c2[3] >> 8 for i in range(3)]
602810
def _call(self, x):
        """"""
        if self.impl == 'pywt':
            coeffs = pywt.wavedecn(
                x, wavelet=self.pywt_wavelet, level=self.nlevels,
                mode=self.pywt_pad_mode, axes=self.axes)
            return pywt.ravel_coeffs(coeffs, axes=self.axes)[0]
        else:
            raise RuntimeError("bad `impl` '{}'".format(self.impl))
378130
def pivot_wavelength_ee(bpass):
    
    from scipy.integrate import simps
    return np.sqrt(simps(bpass.resp, bpass.wlen) /
                   simps(bpass.resp / bpass.wlen**2, bpass.wlen))
143737
def _indiameter_from_fibers(self):
        
        """
        Np = self.num_pores()
        indiam = np.zeros(Np, dtype=float)
        incen = np.zeros([Np, 3], dtype=float)
        hull_pores = np.unique(self._hull_image)
        (Lx, Ly, Lz) = np.shape(self._hull_image)
        (indx, indy, indz) = np.indices([Lx, Ly, Lz])
        indx = indx.flatten()
        indy = indy.flatten()
        indz = indz.flatten()
        for i, pore in enumerate(hull_pores):
            logger.info("Processing pore: "+str(i)+" of "+str(len(hull_pores)))
            dt_pore = self._dt_image*(self._hull_image == pore)
            indiam[pore] = dt_pore.max()*2
            max_ind = np.argmax(dt_pore)
            incen[pore, 0] = indx[max_ind]
            incen[pore, 1] = indy[max_ind]
            incen[pore, 2] = indz[max_ind]
        indiam *= self.network.resolution
        incen *= self.network.resolution
        return (indiam, incen)
1208127
def _make_tonnetz_matrix():
    """"""
    pi = np.pi
    chroma = np.arange(12)

    # Define each row of the transform matrix
    fifth_x = r_fifth*(np.sin((7*pi/6) * chroma))
    fifth_y = r_fifth*(np.cos((7*pi/6) * chroma))
    minor_third_x = r_minor_thirds*(np.sin(3*pi/2 * chroma))
    minor_third_y = r_minor_thirds*(np.cos(3*pi/2 * chroma))
    major_third_x = r_major_thirds*(np.sin(2*pi/3 * chroma))
    major_third_y = r_major_thirds*(np.cos(2*pi/3 * chroma))

    # Return the tonnetz matrix
    return np.vstack((fifth_x, fifth_y,
                      minor_third_x, minor_third_y,
                      major_third_x, major_third_y))
1051311
def formatLij(Lij0,Ne):
	""""""
	#We create Lij as a matrix of lists of laser indices
	global Lij
	Lij=[]
	for i in range(Ne):
		fila=[]
		for j in range(Ne):
			band=False
			for tri in Lij0:
				if [i+1,j+1] == tri[:2]:
					band=True
					break
				elif [j+1,i+1] == tri[:2]:
					band=True
					break
			if band:
				fila+=[tri[2]]
			else:
				fila+=[[]]
		Lij+=[fila]
	return Lij
372213
def get_centrals(self):
        
        result = {}
        for chname in self.rsr.keys():
            result[chname] = {}
            if self.wavespace == "wavelength":
                x__ = self.rsr[chname]["wavelength"]
            else:
                x__ = self.rsr[chname]["wavenumber"]
            for sat in self.rsr[chname].keys():
                if sat in ["wavelength", "wavenumber"]:
                    continue
                if type(self.rsr[chname][sat]) is dict:
                    result[chname][sat] = {}
                    for name in self.rsr[chname][sat].keys():
                        resp = self.rsr[chname][sat][name]
                        result[chname][sat][name] = get_central_wave(x__, resp)
                else:
                    resp = self.rsr[chname][sat]
                    result[chname][sat] = get_central_wave(x__, resp)

        if self.wavespace == "wavelength":
            self.central_wavelength = result
        else:
            self.central_wavenumber = result
1001103
def _complex_even(f, fx, x, h):
        
        n = len(x)
        ee = np.diag(h)
        hess = 2. * np.outer(h, h)
        for i in range(n):
            for j in range(i, n):
                hess[i, j] = (f(x + 1j * ee[i] + ee[j])
                              - f(x + 1j * ee[i] - ee[j])).imag / hess[j, i]
                hess[j, i] = hess[i, j]
        return hess