Gauss gradient

 @param [Float] sigma Spread of Gauss gradient.
 @param [Integer] direction Orientation of Gauss gradient.
 @param [Float] max_error Error of approximated filter.

 @return [Node] Result of filter operation.
435777
def gauss_gradient( sigma, direction, max_error = 1.0 / 0x100 )
      filter_type = DFLOAT.align typecode
      gradient =
        Sequence[ *Array.gauss_gradient_filter( sigma, max_error / dimension ) ].
        to_type filter_type
      blur =
        Sequence[ *Array.gauss_blur_filter( sigma, max_error / dimension ) ].
        to_type filter_type
      ( dimension - 1 ).downto( 0 ).inject self do |retval,i|
        filter = i == direction ? gradient : blur
        retval.convolve filter
      end.force
    end
1263071
def _make_linear_gradient(width, height, angle=90.):
    """"""
    import numpy as np
    X, Y = np.meshgrid(np.linspace(0, 1, width), np.linspace(0, 1, height))
    theta = np.radians(angle % 360)
    c, s = np.cos(theta), np.sin(theta)
    if 0 <= theta and theta < 0.5 * np.pi:
        Z = np.abs(c * X + s * Y)
    elif 0.5 * np.pi <= theta and theta < np.pi:
        Z = np.abs(c * (X - width) + s * Y)
    elif np.pi <= theta and theta < 1.5 * np.pi:
        Z = np.abs(c * (X - width) + s * (Y - height))
    elif 1.5 * np.pi <= theta and theta < 2.0 * np.pi:
        Z = np.abs(c * X + s * (Y - height))
    return (Z - Z.min()) / (Z.max() - Z.min())
47445
def _get_vqsr_training(filter_type, vrn_files, gatk_type):
    
    params = []
    for name, train_info, fname in _get_training_data(vrn_files)[filter_type]:
        if gatk_type == "gatk4":
            params.extend(["--resource:%s,%s" % (name, train_info), fname])
            if filter_type == "INDEL":
                params.extend(["--max-gaussians", "4"])
        else:
            params.extend(["-resource:%s,VCF,%s" % (name, train_info), fname])
            if filter_type == "INDEL":
                params.extend(["--maxGaussians", "4"])
    return params
1038194
def girdle_error(ax, fit, **kwargs):
    
    vertices = []
    codes = []
    for sheet in ('upper','lower'):
        err = plane_errors(fit.axes, fit.covariance_matrix,
                sheet=sheet)
        lonlat = N.array(err)
        lonlat *= -1
        n = len(lonlat)
        if sheet == 'lower':
            lonlat = lonlat[::-1]
        vertices += list(lonlat)
        codes.append(Path.MOVETO)
        codes += [Path.LINETO]*(n-1)

    plot_patch(ax, vertices, codes, **kwargs)
978014
private void computeGradients( float kernelRadius, int kernelWidth ) {

        // generate the gaussian convolution masks
        float kernel[] = new float[kernelWidth];
        float diffKernel[] = new float[kernelWidth];
        int kwidth;
        for( kwidth = 0; kwidth < kernelWidth; kwidth++ ) {
            float g1 = gaussian(kwidth, kernelRadius);
            if (g1 <= GAUSSIAN_CUT_OFF && kwidth >= 2)
                break;
            float g2 = gaussian(kwidth - 0.5f, kernelRadius);
            float g3 = gaussian(kwidth + 0.5f, kernelRadius);
            kernel[kwidth] = (g1 + g2 + g3) / 3f
                    / (2f * (float) Math.PI * kernelRadius * kernelRadius);
            diffKernel[kwidth] = g3 - g2;
        }

        int initX = kwidth - 1;
        int maxX = width - (kwidth - 1);
        int initY = width * (kwidth - 1);
        int maxY = width * (height - (kwidth - 1));

        // perform convolution in x and y directions
        for( int x = initX; x < maxX; x++ ) {
            for( int y = initY; y < maxY; y += width ) {
                int index = x + y;
                float sumX = data[index] * kernel[0];
                float sumY = sumX;
                int xOffset = 1;
                int yOffset = width;
                for( ; xOffset < kwidth; ) {
                    sumY += kernel[xOffset] * (data[index - yOffset] + data[index + yOffset]);
                    sumX += kernel[xOffset] * (data[index - xOffset] + data[index + xOffset]);
                    yOffset += width;
                    xOffset++;
                }

                yConv[index] = sumY;
                xConv[index] = sumX;
            }

        }

        for( int x = initX; x < maxX; x++ ) {
            for( int y = initY; y < maxY; y += width ) {
                float sum = 0f;
                int index = x + y;
                for( int i = 1; i < kwidth; i++ )
                    sum += diffKernel[i] * (yConv[index - i] - yConv[index + i]);

                xGradient[index] = sum;
            }

        }

        for( int x = kwidth; x < width - kwidth; x++ ) {
            for( int y = initY; y < maxY; y += width ) {
                float sum = 0.0f;
                int index = x + y;
                int yOffset = width;
                for( int i = 1; i < kwidth; i++ ) {
                    sum += diffKernel[i] * (xConv[index - yOffset] - xConv[index + yOffset]);
                    yOffset += width;
                }

                yGradient[index] = sum;
            }

        }

        initX = kwidth;
        maxX = width - kwidth;
        initY = width * kwidth;
        maxY = width * (height - kwidth);
        for( int x = initX; x < maxX; x++ ) {
            for( int y = initY; y < maxY; y += width ) {
                int index = x + y;
                int indexN = index - width;
                int indexS = index + width;
                int indexW = index - 1;
                int indexE = index + 1;
                int indexNW = indexN - 1;
                int indexNE = indexN + 1;
                int indexSW = indexS - 1;
                int indexSE = indexS + 1;

                float xGrad = xGradient[index];
                float yGrad = yGradient[index];
                float gradMag = hypot(xGrad, yGrad);

                // perform non-maximal supression
                float nMag = hypot(xGradient[indexN], yGradient[indexN]);
                float sMag = hypot(xGradient[indexS], yGradient[indexS]);
                float wMag = hypot(xGradient[indexW], yGradient[indexW]);
                float eMag = hypot(xGradient[indexE], yGradient[indexE]);
                float neMag = hypot(xGradient[indexNE], yGradient[indexNE]);
                float seMag = hypot(xGradient[indexSE], yGradient[indexSE]);
                float swMag = hypot(xGradient[indexSW], yGradient[indexSW]);
                float nwMag = hypot(xGradient[indexNW], yGradient[indexNW]);
                float tmp;
                /*
                 * An explanation of what's happening here, for those who want
                 * to understand the source: This performs the "non-maximal
                 * supression" phase of the Canny edge detection in which we
                 * need to compare the gradient magnitude to that in the
                 * direction of the gradient; only if the value is a local
                 * maximum do we consider the point as an edge candidate.
                 * 
                 * We need to break the comparison into a number of different
                 * cases depending on the gradient direction so that the
                 * appropriate values can be used. To avoid computing the
                 * gradient direction, we use two simple comparisons: first we
                 * check that the partial derivatives have the same sign (1)
                 * and then we check which is larger (2). As a consequence, we
                 * have reduced the problem to one of four identical cases that
                 * each test the central gradient magnitude against the values at
                 * two points with 'identical support'; what this means is that
                 * the geometry required to accurately interpolate the magnitude
                 * of gradient function at those points has an identical
                 * geometry (upto right-angled-rotation/reflection).
                 * 
                 * When comparing the central gradient to the two interpolated
                 * values, we avoid performing any divisions by multiplying both
                 * sides of each inequality by the greater of the two partial
                 * derivatives. The common comparand is stored in a temporary
                 * variable (3) and reused in the mirror case (4).
                 * 
                 */
                if (xGrad * yGrad <= (float) 0 /*(1)*/
                ? Math.abs(xGrad) >= Math.abs(yGrad) /*(2)*/
                ? (tmp = Math.abs(xGrad * gradMag)) >= Math.abs(yGrad * neMag - (xGrad + yGrad)
                        * eMag) /*(3)*/
                        && tmp > Math.abs(yGrad * swMag - (xGrad + yGrad) * wMag) /*(4)*/
                : (tmp = Math.abs(yGrad * gradMag)) >= Math.abs(xGrad * neMag - (yGrad + xGrad)
                        * nMag) /*(3)*/
                        && tmp > Math.abs(xGrad * swMag - (yGrad + xGrad) * sMag) /*(4)*/
                : Math.abs(xGrad) >= Math.abs(yGrad) /*(2)*/
                ? (tmp = Math.abs(xGrad * gradMag)) >= Math.abs(yGrad * seMag + (xGrad - yGrad)
                        * eMag) /*(3)*/
                        && tmp > Math.abs(yGrad * nwMag + (xGrad - yGrad) * wMag) /*(4)*/
                : (tmp = Math.abs(yGrad * gradMag)) >= Math.abs(xGrad * seMag + (yGrad - xGrad)
                        * sMag) /*(3)*/
                        && tmp > Math.abs(xGrad * nwMag + (yGrad - xGrad) * nMag) /*(4)*/
                ) {
                    magnitude[index] = gradMag >= MAGNITUDE_LIMIT
                            ? MAGNITUDE_MAX
                            : (int) (MAGNITUDE_SCALE * gradMag);
                    // NOTE: The orientation of the edge is not employed by this
                    // implementation. It is a simple matter to compute it at
                    // this point as: Math.atan2(yGrad, xGrad);
                } else {
                    magnitude[index] = 0;
                }
            }
        }
    }
1232612
def adjust_distributions_using_gradient
      return if @recent_scores.length < 3
      vx = (1..@recent_scores.length).to_a.to_numeric
      vy = @recent_scores.reverse.to_numeric
      r = Statsample::Regression::Simple.new_from_vectors(vx,vy)
      slope = r.b
      if slope > 0
        @distributions.each_pair { |k, d| d.tighten slope }
      elsif slope < 0
        @distributions.each_pair { |k, d| d.loosen slope }
      end
    end
110284
private double computeOutputError(double output, double[] gradient) {

        double error = 0.0;
        double out = outputLayer.output[0];
        double g = output - out;

        error += (0.5*g * g);


        gradient[0] = g;

        return error;
    }
1321824
def _validI(x, y, weights):
    '''
     
    '''
    # density filter:
    i = np.logical_and(np.isfinite(y), weights > np.median(weights))
    # filter outliers:
    try:
        grad = np.abs(np.gradient(y[i]))
        max_gradient = 4 * np.median(grad)
        i[i][grad > max_gradient] = False
    except (IndexError, ValueError):
        pass
    return i
1021625
def gauss_fltr_opencv(dem, size=3, sigma=1):
    
    import cv2
    dem = malib.checkma(dem)
    dem_cv = cv2.GaussianBlur(dem.filled(np.nan), (size, size), sigma)
    out = np.ma.fix_invalid(dem_cv)
    out.set_fill_value(dem.fill_value)
    return out
305269
def get_smeared_densities(self, sigma):
        

        from scipy.ndimage.filters import gaussian_filter1d
        diff = [self.frequencies[i + 1] - self.frequencies[i]
                for i in range(len(self.frequencies) - 1)]
        avgdiff = sum(diff) / len(diff)

        smeared_dens = gaussian_filter1d(self.densities, sigma / avgdiff)
        return smeared_dens
588444
def gaussian_smooth_eval(*vals, sigma_vec)
		npix = 7
		raise "npix must be odd" if npix%2==0
		case vals.size
		when 2
# 			delt0 = 3.0*0.999999*sigma_vec[0] / (npix-1)
# 			delt1 = 3.0*0.999999*sigma_vec[1] / (npix-1)
# 			sig3 = 3.0*sigma
			vals0 = GSL::Vector.linspace(vals[0] - 3.0* sigma_vec[0], vals[0] + 3.0* sigma_vec[0], npix)
			vals1 = GSL::Vector.linspace(vals[1] - 3.0* sigma_vec[1], vals[1] + 3.0* sigma_vec[1], npix)
			mat = GSL::Matrix.alloc(vals0.size, vals1.size)
			for i in 0...vals0.size
				for j in 0...vals1.size
					mat[i,j] = eval(vals0[i], vals1[j])
				end
			end
			mat.gaussian_smooth(*sigma_vec.to_a)
			cent = (npix - 1) / 2
			return mat[cent, cent]
			
		else
			raise 'Not supported for this number of dimensions yet'
		end
	end
37811
protected double computeErrors(Se3_F64 fiducialToCamera ) {
		if( fiducialToCamera.T.z < 0 ) {
			// the low level algorithm should already filter this code, but just incase
			return Double.MAX_VALUE;
		}

		double maxError = 0;

		for( int i = 0; i < 4; i++ ) {
			maxError = Math.max(maxError,computePixelError(fiducialToCamera, points.get(i).location, listObs.get(i)));
		}

		return maxError;
	}
913529
def make_gaussian_kernel(sigma, npix=501, cdelt=0.01, xpix=None, ypix=None):
    

    sigma /= cdelt

    def fn(t, s): return 1. / (2 * np.pi * s ** 2) * np.exp(
        -t ** 2 / (s ** 2 * 2.0))
    dxy = make_pixel_distance(npix, xpix, ypix)
    k = fn(dxy, sigma)
    k /= (np.sum(k) * np.radians(cdelt) ** 2)

    return k
325333
public static <T extends ImageBase<T>> BlurStorageFilter<T> gaussian(ImageType<T> type , double sigma , int radius ) {
		return new BlurStorageFilter<>("gaussian", type, sigma, radius);
	}
1222534
def vel_disp(self, kwargs_profile, kwargs_aperture, kwargs_light, kwargs_anisotropy, num=1000):
        
        sigma_s2_sum = 0
        for i in range(0, num):
            sigma_s2_draw = self._vel_disp_one(kwargs_profile, kwargs_aperture, kwargs_light, kwargs_anisotropy)
            sigma_s2_sum += sigma_s2_draw
        sigma_s2_average = sigma_s2_sum/num
        return np.sqrt(sigma_s2_average)
855827
def _rotation_direction_in_detector(geometry):
    """"""
    du, dv = geometry.det_axes_init
    axis = geometry.axis
    det_normal = np.cross(dv, du)
    rot_dir = np.cross(axis, det_normal)
    c = np.array([np.vdot(rot_dir, du), np.vdot(rot_dir, dv)])
    cnorm = np.linalg.norm(c)

    # Check for numerical errors
    assert cnorm != 0

    return c / cnorm
1316092
def bound_perturbed_gmres(pseudo, p, epsilon, deltas):
    '''
    '''
    if not numpy.all(numpy.array(deltas) > epsilon):
        raise ArgumentError('all deltas have to be greater than epsilon')

    bound = []
    for delta in deltas:
        # get boundary paths
        paths = pseudo.contour_paths(delta)

        # get vertices on boundary
        vertices = paths.vertices()

        # evaluate polynomial
        supremum = numpy.max(numpy.abs(p(vertices)))

        # compute bound
        bound.append(epsilon/(delta-epsilon)
                     * paths.length()/(2*numpy.pi*delta)
                     * supremum)
    return bound