我们从Python开源项目中，提取了以下11个代码示例，用于说明如何使用scipy.ndimage.mean()。
def apply(self): sc = self.scene org = self.original factor = self.factor sx, sy = sc.displacement.shape gx, gy = num.ogrid[0:sx, 0:sy] regions = sy/factor * (gx/factor) + gy/factor indices = num.arange(regions.max() + 1) def block_downsample(arr): res = ndimage.mean( arr, labels=regions, index=indices) res.shape = (sx/factor, sy/factor) return res sc.displacement = block_downsample(sc.displacement) sc.theta = block_downsample(sc.theta) sc.phi = block_downsample(sc.phi) sc.frame.dLat = org['frame.dLat'] * self.factor sc.frame.dLon = org['frame.dLat'] * self.factor
def downsample(voxels, step, method='max'): """ downsample a voxels matrix by a factor of step. downsample method options: max/mean same as a pooling """ assert step > 0 assert voxels.ndim == 3 or voxels.ndim == 4 assert method in ('max', 'mean') if step == 1: return voxels if voxels.ndim == 3: sx, sy, sz = voxels.shape[-3:] X, Y, Z = np.ogrid[0:sx, 0:sy, 0:sz] regions = sz/step * sy/step * (X/step) + sz/step * (Y/step) + Z/step if method == 'max': res = ndimage.maximum(voxels, labels=regions, index=np.arange(regions.max() + 1)) elif method == 'mean': res = ndimage.mean(voxels, labels=regions, index=np.arange(regions.max() + 1)) res.shape = (sx/step, sy/step, sz/step) return res else: res0 = downsample(voxels[0], step, method) res = np.zeros((voxels.shape[0],) + res0.shape) res[0] = res0 for ind in xrange(1, voxels.shape[0]): res[ind] = downsample(voxels[ind], step, method) return res
def property_topomesh_area_smoothing_force(topomesh,target_areas=None): compute_topomesh_property(topomesh,'vertices',degree=2) compute_topomesh_property(topomesh,'length',degree=1) triangle_vertices = topomesh.wisp_property('vertices',degree=2).values() rotated_triangle_vertices = np.transpose([triangle_vertices[:,2],triangle_vertices[:,0],triangle_vertices[:,1]]) antirotated_triangle_vertices = np.transpose([triangle_vertices[:,1],triangle_vertices[:,2],triangle_vertices[:,0]]) triangle_vertices = np.append(np.append(triangle_vertices,rotated_triangle_vertices,axis=0),antirotated_triangle_vertices,axis=0) edge_index_list = np.array([[1, 2],[0, 1],[0, 2]]) triangle_edge_vertices = triangle_vertices[:,edge_index_list] triangle_edge_vectors = topomesh.wisp_property('barycenter',degree=0).values(triangle_edge_vertices[...,1]) - topomesh.wisp_property('barycenter',degree=0).values(triangle_edge_vertices[...,0]) #triangle_edge_lengths = np.power(np.sum(np.power(triangle_edge_vectors,2.0),axis=2),0.5) triangle_edge_lengths = np.linalg.norm(triangle_edge_vectors,axis=2) triangle_edge_directions = triangle_edge_vectors/triangle_edge_lengths[...,np.newaxis] triangle_perimeters = np.sum(triangle_edge_lengths,axis=1) triangle_areas = np.sqrt((triangle_perimeters/2.0)*(triangle_perimeters/2.0-triangle_edge_lengths[:,0])*(triangle_perimeters/2.0-triangle_edge_lengths[:,1])*(triangle_perimeters/2.0-triangle_edge_lengths[:,2])) triangle_areas[np.where(triangle_areas==0.0)] = 0.001 if target_areas is None: target_areas = triangle_areas.mean()*np.ones(len(triangle_areas)) else: target_areas = np.tile(target_areas,(3)) area_edge_1_force = (target_areas-triangle_areas)*(triangle_edge_lengths[:,0]**2+triangle_edge_lengths[:,2]**2-triangle_edge_lengths[:,1]**2)*triangle_edge_lengths[:,1] area_edge_2_force = (target_areas-triangle_areas)*(triangle_edge_lengths[:,0]**2+triangle_edge_lengths[:,1]**2-triangle_edge_lengths[:,2]**2)*triangle_edge_lengths[:,2] area_unitary_force = -(area_edge_1_force[:,np.newaxis]*triangle_edge_directions[:,1] + area_edge_2_force[:,np.newaxis]*triangle_edge_directions[:,2]) triangle_force = np.transpose([nd.sum(area_unitary_force[:,d],triangle_vertices[:,0],index=list(topomesh.wisps(0))) for d in xrange(3)]) return triangle_force
def run(self, ips, snap, img, para = None): intenimg = WindowsManager.get(para['inten']).ips.img strc = ndimage.generate_binary_structure(2, 1 if para['con']=='4-connect' else 2) buf, n = ndimage.label(snap, strc, output=np.uint16) index = range(1, n+1) idx = (np.ones(n+1)*para['front']).astype(np.uint8) msk = np.ones(n, dtype=np.bool) if para['mean']>0: msk *= ndimage.mean(intenimg, buf, index)>=para['mean'] if para['mean']<0: msk *= ndimage.mean(intenimg, buf, index)<-para['mean'] if para['max']>0: msk *= ndimage.maximum(intenimg, buf, index)>=para['max'] if para['max']<0: msk *= ndimage.maximum(intenimg, buf, index)<-para['max'] if para['min']>0: msk *= ndimage.minimum(intenimg, buf, index)>=para['min'] if para['min']<0: msk *= ndimage.minimum(intenimg, buf, index)<-para['min'] if para['sum']>0: msk *= ndimage.sum(intenimg, buf, index)>=para['sum'] if para['sum']<0: msk *= ndimage.sum(intenimg, buf, index)<-para['sum'] if para['std']>0: msk *= ndimage.standard_deviation(intenimg, buf, index)>=para['std'] if para['std']<0: msk *= ndimage.standard_deviation(intenimg, buf, index)<-para['std'] xy = ndimage.center_of_mass(intenimg, buf, index) xy = np.array(xy).round(2).T idx[1:][~msk] = para['back'] idx[0] = 0 img[:] = idx[buf] WindowsManager.get(para['inten']).ips.mark = RGMark((xy.T, msk)) WindowsManager.get(para['inten']).ips.update = True
def _roundness(x, y, pixelSize): """ Implements equation (5) of Cappellari & Copin (2003) """ n = x.size equivalentRadius = np.sqrt(n/np.pi)*pixelSize xBar, yBar = np.mean(x), np.mean(y) # Geometric centroid here! maxDistance = np.sqrt(np.max((x - xBar)**2 + (y - yBar)**2)) roundness = maxDistance/equivalentRadius - 1. return roundness #----------------------------------------------------------------------
def _reassign_bad_bins(classe, x, y): """ Implements steps (vi)-(vii) in section 5.1 of Cappellari & Copin (2003) """ # Find the centroid of all successful bins. # CLASS = 0 are unbinned pixels which are excluded. # good = np.unique(classe[classe > 0]) xnode = ndimage.mean(x, labels=classe, index=good) ynode = ndimage.mean(y, labels=classe, index=good) # Reassign pixels of bins with S/N < targetSN # to the closest centroid of a good bin # bad = classe == 0 index = voronoi_tessellation(x[bad], y[bad], xnode, ynode, [1]) classe[bad] = good[index] # Recompute all centroids of the reassigned bins. # These will be used as starting points for the CVT. # good = np.unique(classe) xnode = ndimage.mean(x, labels=classe, index=good) ynode = ndimage.mean(y, labels=classe, index=good) return xnode, ynode #----------------------------------------------------------------------------
def property_topomesh_laplacian_epidermis_convexity_force(topomesh): """todo""" if not topomesh.has_wisp_property('vertices',degree=1,is_computed=True): compute_topomesh_property(topomesh,'vertices',degree=1) if not topomesh.has_wisp_property('barycenter',degree=3,is_computed=True): compute_topomesh_property(topomesh,'barycenter',degree=3) if not topomesh.has_wisp_property('epidermis',degree=1,is_computed=True): compute_topomesh_property(topomesh,'epidermis',degree=1) if not topomesh.has_wisp_property('cells',degree=0,is_computed=True): compute_topomesh_property(topomesh,'cells',degree=0) if not topomesh.has_wisp_property('cells',degree=1,is_computed=True): compute_topomesh_property(topomesh,'cells',degree=1) epidermis_convexity_force = array_dict(np.zeros_like(topomesh.wisp_property('barycenter',degree=0).values(),np.float32),np.array(list(topomesh.wisps(0)))) edge_cells_degree = np.array(map(len,topomesh.wisp_property('cells',degree=1).values())) # edge_vertices = topomesh.wisp_property('vertices',degree=1).values()[np.where((topomesh.wisp_property('epidermis',degree=1).values()))] # edge_vertices = topomesh.wisp_property('vertices',degree=1).values()[np.where((topomesh.wisp_property('epidermis',degree=1).values())&(edge_cells_degree>1))] edge_vertices = topomesh.wisp_property('vertices',degree=1).values()[np.where((edge_cells_degree>2)|((topomesh.wisp_property('epidermis',degree=1).values())&(edge_cells_degree>1)))] epidermis_vertices = np.unique(edge_vertices) vertices_degrees = np.array(nd.sum(np.ones_like(edge_vertices),edge_vertices,index=epidermis_vertices),int) reversed_edge_vertices = np.transpose([edge_vertices[:,1],edge_vertices[:,0]]) edge_vertices = np.append(edge_vertices,reversed_edge_vertices,axis=0) edge_vectors = topomesh.wisp_property('barycenter',degree=0).values(edge_vertices[:,1]) - topomesh.wisp_property('barycenter',degree=0).values(edge_vertices[:,0]) laplacian_force = np.transpose([nd.sum(edge_vectors[:,0],edge_vertices[:,0],index=epidermis_vertices), nd.sum(edge_vectors[:,1],edge_vertices[:,0],index=epidermis_vertices), nd.sum(edge_vectors[:,2],edge_vertices[:,0],index=epidermis_vertices)]) laplacian_force = laplacian_force/vertices_degrees[:,np.newaxis] # vertex_cell_barycenter = np.array([np.mean(topomesh.wisp_property('barycenter',degree=3).values(c),axis=0) for c in topomesh.wisp_property('cells',0).values(epidermis_vertices)]) # vertices_directions = topomesh.wisp_property('barycenter',degree=0).values(epidermis_vertices) - vertex_cell_barycenter # vertices_directions = vertices_directions/np.linalg.norm(vertices_directions,axis=1)[:,np.newaxis] # vertices_products = np.einsum('ij,ij->i',laplacian_force,vertices_directions) # convex_points = np.where(vertices_products<0.)[0] # laplacian_force[convex_points] -= vertices_directions[convex_points]*vertices_products[convex_points,np.newaxis] # laplacian_force[convex_points] = np.array([0.,0.,0.]) epidermis_convexity_force.update(laplacian_force,keys=epidermis_vertices,erase_missing_keys=False) return epidermis_convexity_force.values()
def run(self, ips, imgs, para = None): inten = WindowsManager.get(para['inten']).ips if not para['slice']: imgs = [inten.img] msks = [ips.img] else: msks = ips.imgs if len(msks)==1: msks *= len(imgs) buf = imgs[0].astype(np.uint16) strc = ndimage.generate_binary_structure(2, 1 if para['con']=='4-connect' else 2) idct = ['Max','Min','Mean','Variance','Standard','Sum'] key = {'Max':'max','Min':'min','Mean':'mean', 'Variance':'var','Standard':'std','Sum':'sum'} idct = [i for i in idct if para[key[i]]] titles = ['Slice', 'ID'][0 if para['slice'] else 1:] if para['center']: titles.extend(['Center-X','Center-Y']) if para['extent']: titles.extend(['Min-Y','Min-X','Max-Y','Max-X']) titles.extend(idct) k = ips.unit[0] data, mark = [], [] for i in range(len(imgs)): n = ndimage.label(msks[i], strc, output=buf) index = range(1, n+1) dt = [] if para['slice']:dt.append([i]*n) dt.append(range(n)) xy = ndimage.center_of_mass(imgs[i], buf, index) xy = np.array(xy).round(2).T if para['center']:dt.extend([xy[1]*k, xy[0]*k]) boxs = [None] * n if para['extent']: boxs = ndimage.find_objects(buf) boxs = [(i[0].start, i[1].start, i[0].stop, i[1].stop) for i in boxs] for j in (0,1,2,3): dt.append([i[j]*k for i in boxs]) if para['max']:dt.append(ndimage.maximum(imgs[i], buf, index).round(2)) if para['min']:dt.append(ndimage.minimum(imgs[i], buf, index).round(2)) if para['mean']:dt.append(ndimage.mean(imgs[i], buf, index).round(2)) if para['var']:dt.append(ndimage.variance(imgs[i], buf, index).round(2)) if para['std']:dt.append(ndimage.standard_deviation(imgs[i], buf, index).round(2)) if para['sum']:dt.append(ndimage.sum(imgs[i], buf, index).round(2)) mark.append([(center, cov) for center,cov in zip(xy.T, boxs)]) data.extend(list(zip(*dt))) IPy.table(inten.title+'-region statistic', data, titles) inten.mark = Mark(mark) inten.update = True
def _cvt_equal_mass(x, y, signal, noise, xnode, ynode, pixelsize, quiet, sn_func, wvt): """ Implements the modified Lloyd algorithm in section 4.1 of Cappellari & Copin (2003). NB: When the keyword WVT is set this routine includes the modification proposed by Diehl & Statler (2006). """ dens2 = (signal/noise)**4 # See beginning of section 4.1 of CC03 scale = np.ones_like(xnode) # Start with the same scale length for all bins for it in range(1, xnode.size): # Do at most xnode.size iterations xnode_old, ynode_old = xnode.copy(), ynode.copy() classe = voronoi_tessellation(x, y, xnode, ynode, scale) # Computes centroids of the bins, weighted by dens**2. # Exponent 2 on the density produces equal-mass Voronoi bins. # The geometric centroids are computed if WVT keyword is set. # good = np.unique(classe) if wvt: for k in good: index = np.flatnonzero(classe == k) # Find subscripts of pixels in bin k. xnode[k], ynode[k] = np.mean(x[index]), np.mean(y[index]) sn = sn_func(index, signal, noise) scale[k] = np.sqrt(index.size/sn) # Eq. (4) of Diehl & Statler (2006) else: mass = ndimage.sum(dens2, labels=classe, index=good) xnode = ndimage.sum(x*dens2, labels=classe, index=good)/mass ynode = ndimage.sum(y*dens2, labels=classe, index=good)/mass diff2 = np.sum((xnode - xnode_old)**2 + (ynode - ynode_old)**2) diff = np.sqrt(diff2)/pixelsize if not quiet: print('Iter: %4i Diff: %.4g' % (it, diff)) if diff < 0.1: break # If coordinates have changed, re-compute (Weighted) Voronoi Tessellation of the pixels grid # if diff > 0: classe = voronoi_tessellation(x, y, xnode, ynode, scale) good = np.unique(classe) # Check for zero-size Voronoi bins # Only return the generators and scales of the nonzero Voronoi bins return xnode[good], ynode[good], scale[good], it #-----------------------------------------------------------------------
