importnumpyasnpimportmatplotlib.pyplotaspltimportcv2fromcollectionsimportdefaultdictfromscipyimportndimageclassChannel():colormaps={"k":"Greys","r":"Reds","g":"Greens","b":"Blues"}def__init__(self,data=None):ifdataisNone:self.data=np.zeros((0,0),dtype=np.uint8)else:self.data=self.import_array(data)def__call__(self,data=None):""""Set and return current channel data."""ifnot(dataisNone):self.data=self.import_array(data)returnself.datadefimport_array(self,data):data=np.where(data<0,0,data)data=np.where(data>255,255,data)returndata.astype(np.uint8)@propertydefhistogram(self):returncv2.calcHist([self.data],[0],None,[256],[0,255]).flatten()defshow(self,ax,color='k'):ax.imshow(self.data,cmap=self.colormaps[color])defshow_histogram(self,ax,color='k'):bins=np.arange(0,256)ax.fill_between(bins,self.histogram,0,color=color,alpha=0.3)defequalize_histogram(self):self.data=cv2.equalizeHist(self.data)defblur_Gaussian(self,ksize=3):self.data=cv2.GaussianBlur(self.data,(ksize,ksize),0)defthreshold(self,T):T,self.data=cv2.threshold(self.data,T,255,cv2.THRESH_BINARY)defthreshold_otsu(self):T,self.data=cv2.threshold(self.data,0,255,cv2.THRESH_BINARY+cv2.THRESH_OTSU)returnTdefmorph_close(self,iterations=1,kernel=np.ones((3,3))):ret=cv2.dilate(self.data,kernel=kernel,iterations=iterations)self.data=cv2.erode(ret,kernel=kernel,iterations=iterations)defmorph_open(self,iterations=1,kernel=np.ones((3,3))):ret=cv2.erode(self.data,kernel=kernel,iterations=iterations)self.data=cv2.dilate(ret,kernel=kernel,iterations=iterations)deffind_particles(self):structure=np.ones((3,3))labels,particle_count=ndimage.measurements.label(self.data,structure=structure)# locate centers of masscenters_of_mass=ndimage.center_of_mass(self.data,labels,np.arange(1,particle_count+1))# get slices containting particlesslices=ndimage.find_objects(labels)# get raster points for each particlepoints=[]Y,X=labels.shapeforninrange(len(centers_of_mass)):locs=labels==(n+1)points.append(np.array([(x,y)forxinrange(X)foryinrange(Y)iflocs[y,x]]))returncenters_of_mass,slices,pointsclassParticle():def__init__(self,center_of_mass,slice,raster_points):self.yc,self.xc=center_of_massyslice,xslice=sliceself.x1,self.y1=xslice.start,yslice.startself.x2,self.y2=xslice.stop,yslice.stopself.w=self.x2-self.x1self.h=self.y2-self.y1self.size=np.sqrt(self.w*self.h)self.raster_points=raster_pointsdefshow(self,ax):ya=self.y1yb=self.y2xa=self.x1xb=self.x2ax.plot([xa,xb,xb,xa,xa],[yb,yb,ya,ya,yb],'r',lw=1)classLabeler():def__init__(self):self.channels=defaultdict(Channel)self.particles=[]defread(self,filepath):self.filepath=filepathself.img=cv2.cvtColor(cv2.imread(filepath),cv2.COLOR_BGR2RGB)defshow(self,ax):ax.imshow(self.img)ax.set_title(self.filepath)ax.set_xlabel("<--- w --->")ax.set_ylabel("<--- h --->")forpinlabeler.particles:p.show(ax)deffind_particles(self,filepath,size_range=(0,1e6)):self.read(filepath)r=self.img[:,:,0]g=self.img[:,:,1]b=self.img[:,:,2]self.channels['bw'](g-0.5*b-0.1*r-20.0)self.channels['bw'].equalize_histogramself.channels['bw'].blur_Gaussian(31)self.channels['bw'].threshold_otsu()self.channels['bw'].morph_close(3)self.channels['bw'].morph_open(8)centers_of_mass,slices,raster_points=self.channels['bw'].find_particles()self.particles=[Particle(centers_of_mass[k],slices[k],raster_points[k])forkinrange(len(centers_of_mass))]self.particles=[pforpinself.particlesif(p.size>=size_range[0])and(p.size<=size_range[1])]