analytic solution of integral for NFW profile to compute deflection angel and gamma

        :param x: R/Rs
        :type x: float >0
343163
def _g(self, x, tau):
        
        return tau ** 2 * (tau ** 2 + 1) ** -2 * (
                (tau ** 2 + 1 + 2 * (x ** 2 - 1)) * self.F(x) + tau * np.pi + (tau ** 2 - 1) * np.log(tau) +
                np.sqrt(tau ** 2 + x ** 2) * (-np.pi + self.L(x, tau) * (tau ** 2 - 1) * tau ** -1))
1042188
def cnst_A(self, X, Xf=None):
        
        """

        if Xf is None:
            Xf = sl.rfftn(X, axes=self.axes)
        return sl.irfftn(self.GAf*Xf[..., np.newaxis], self.axsz,
                         axes=self.axes)
843197
def _F(self, x, p):

        
        prefactor = 0.5 * (1 + p ** 2) ** -1 * p

        if isinstance(x, np.ndarray):

            inds0 = np.where(x * p == 1)
            inds1 = np.where(x * p < 1)
            inds2 = np.where(x * p > 1)

            func = np.ones_like(x)

            func[inds0] = self._u(x[inds0]) ** -1 * (np.pi + 2 * p * np.arctanh(self._u(x[inds0]) ** -1))

            func[inds1] = self._u(x[inds1]) ** -1 * (np.pi + 2 * p * np.arctanh(self._u(x[inds1]) ** -1)) - \
                          (2 * p * self._g(x[inds1], p) ** -1 * np.arctanh(self._g(x[inds1], p)))

            func[inds2] = self._u(x[inds2]) ** -1 * (np.pi + 2 * p * np.arctanh(self._u(x[inds2]) ** -1)) - \
                          (2 * p * self._f(x[inds2], p) ** -1 * np.arctan(self._f(x[inds2], p)))

            return prefactor * func

        else:

            if x * p == 1:
                func = self._u(x) ** -1 * (np.pi + 2 * p * np.arctanh(self._u(x) ** -1))
            elif x * p < 1:
                func = self._u(x) ** -1 * (np.pi + 2 * p * np.arctanh(self._u(x) ** -1)) - \
                       (2 * p * self._g(x, p) ** -1 * np.arctanh(self._g(x, p)))
            else:
                func = self._u(x) ** -1 * (np.pi + 2 * p * np.arctanh(self._u(x) ** -1)) - \
                       (2 * p * self._f(x, p) ** -1 * np.arctan(self._f(x, p)))

            return prefactor * func
912894
def obfn_dfd(self):
        
        """

        Ef = sl.inner(self.Df, self.obfn_fvarf(), axis=self.cri.axisM) - \
            self.Sf
        return sl.rfl2norm2(Ef, self.S.shape, axis=self.cri.axisN) / 2.0
108202
def obfn_dfd(self):
        
        """

        Ef = sl.inner(self.Zf, self.obfn_fvarf(), axis=self.cri.axisM) \
          - self.Sf
        return (np.linalg.norm(self.W * sl.irfftn(Ef, self.cri.Nv,
                                                  self.cri.axisN))**2) / 2.0
483599
def _hammer_function(self, x,x0,r_x0, s_x0):
        '''
        
        '''
        return norm.logcdf((np.sqrt((np.square(np.atleast_2d(x)[:,None,:]-np.atleast_2d(x0)[None,:,:])).sum(-1))- r_x0)/s_x0)
794461
def rn(x, af, rate):
    
    (taus,devs,errs,ns) = at.adev(x,taus=[af*rate], data_type='phase', rate=rate) 
    oadev_x = devs[0]
    (mtaus,mdevs,errs,ns) = at.mdev(x,taus=[af*rate], data_type='phase', rate=rate)
    mdev_x = mdevs[0]
    rn = pow(mdev_x/oadev_x,2)
    return rn
1181830
def prox_hard(X, step, thresh=0):
    
    thresh_ = _step_gamma(step, thresh)
    below = np.abs(X) < thresh_
    X[below] = 0
    return X
1268201
def _F(self, X, b, c = 0.001):
        

        if b == 1:
            b = 1 + c

        prefac = (b - 1) ** -2

        if isinstance(X, np.ndarray):

            X[np.where(X == 1)] = 1 - c

            output = np.empty_like(X)

            inds1 = np.where(np.absolute(X - b)<c)
            output[inds1] = prefac*(-2 - b + (1 + b + b ** 2) * self._nfw_func(b)) * (1 + b) ** -1

            inds2 = np.where(np.absolute(X - b)>=c)

            output[inds2] = prefac * ((X[inds2] ** 2 - 1) ** -1 * (1 - b -
                                    (1 - b * X[inds2] ** 2) * self._nfw_func(X[inds2])) - \
                                        self._nfw_func(X[inds2] * b ** -1))

        else:

            if X == 1:
                X = 1-c

            if np.absolute(X - b)<c:
                output = prefac * (-2 - b + (1 + b + b ** 2) * self._nfw_func(b)) * (1 + b) ** -1

            else:
                output = prefac * ((X ** 2 - 1) ** -1 * (1 - b -
                                    (1 - b * X ** 2) * self._nfw_func(X)) - \
                                        self._nfw_func(X * b ** -1))

        return output
709417
def _optm(x, alpha, mu, sigma):
    '''
    '''
    return ss.skewnorm.pdf(x, alpha, mu, sigma)
1267211
def xstep(self):
        
        """

        b = self.AHSf + self.rho*np.sum(
            np.conj(self.Gf)*sl.rfftn(self.Y-self.U, axes=self.axes),
            axis=self.Y.ndim-1)
        self.Xf = b / (self.AHAf + self.rho*self.GHGf)
        self.X = sl.irfftn(self.Xf, self.axsz, axes=self.axes)

        if self.opt['LinSolveCheck']:
            ax = (self.AHAf + self.rho*self.GHGf)*self.Xf
            self.xrrs = sl.rrs(ax, b)
        else:
            self.xrrs = None
956023
def _call(self, x):
        """"""

        # Convert to array with most
        arr = self._asarray(x)
        svd_diag = np.linalg.svd(arr, compute_uv=False)

        # Rotate the axes so the svd-direction is first
        s_reordered = moveaxis(svd_diag, -1, 0)

        # Return nuclear norm
        return self.outernorm(self.pwisenorm(s_reordered))
1253441
def calcdelay(self, ant1, ant2, skyfreq, pol):
        

        select = self.select[n.where( (self.skyfreq[self.select] == skyfreq) & (self.polarization[self.select] == pol) )[0]]

        ind1 = n.where(ant1 == self.antnum[select])
        ind2 = n.where(ant2 == self.antnum[select])
        d1 = self.delay[select][ind1]
        d2 = self.delay[select][ind2]
        if len(d1-d2) > 0:
            return d1-d2
        else:
            return n.array([0])
989047
def expected_magx(RAW_IMU, ATTITUDE, inclination, declination):
    ''''''
    v = expected_mag(RAW_IMU, ATTITUDE, inclination, declination)
    return v.x
1374683
def _call(self, x):
        """"""

        sum_constr = abs(x.ufuncs.sum() / self.diameter - 1) <= self.sum_rtol

        nonneq_constr = x.ufuncs.greater_equal(0).asarray().all()

        if sum_constr and nonneq_constr:
            return 0
        else:
            return np.inf
786755
def fit_deriv(x, amplitude, mean, stddev):
        
        import operator
        return list(map(
            operator.neg, Gaussian1D.fit_deriv(x, amplitude, mean, stddev)))
951874
def xAxisIsMajor(self):
        '''
        
        '''
        return max(self.radius.x, self.radius.y) == self.radius.x