Mercurial > repos > perssond > coreograph

comparison UNetCoreograph.py @ 1:57f1260ca94e draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
"planemo upload commit fec9dc76b3dd17b14b02c2f04be9d30f71eba1ae"
author watsocam
date Fri, 11 Mar 2022 23:40:51 +0000
parents 99308601eaa6
children
comparison
equal deleted inserted replaced
0:99308601eaa6 1:57f1260ca94e
1 import numpy as np 1 import numpy as np
2 from scipy import misc as sm 2 from scipy import misc as sm
3 import shutil 3 import shutil
4 import scipy.io as sio 4 import scipy.io as sio
5 import os 5 import os
6 os.environ['TF_CPP_MIN_LOG_LEVEL'] = '2'
7 import logging
8 logging.getLogger('tensorflow').setLevel(logging.FATAL)
6 import skimage.exposure as sk 9 import skimage.exposure as sk
7 import cv2 10 import cv2
8 import argparse 11 import argparse
9 import pytiff 12 import pytiff
10 import tifffile 13 import tifffile
12 from skimage.morphology import * 15 from skimage.morphology import *
13 from skimage.exposure import rescale_intensity 16 from skimage.exposure import rescale_intensity
14 from skimage.segmentation import chan_vese, find_boundaries, morphological_chan_vese 17 from skimage.segmentation import chan_vese, find_boundaries, morphological_chan_vese
15 from skimage.measure import regionprops,label, find_contours 18 from skimage.measure import regionprops,label, find_contours
16 from skimage.transform import resize 19 from skimage.transform import resize
17 from skimage.filters import gaussian 20 from skimage.filters import gaussian, threshold_otsu
18 from skimage.feature import peak_local_max,blob_log 21 from skimage.feature import peak_local_max,blob_log
19 from skimage.color import label2rgb 22 from skimage.color import gray2rgb as gray2rgb
20 import skimage.io as skio 23 import skimage.io as skio
24 from scipy.ndimage.morphology import binary_fill_holes
21 from skimage import img_as_bool 25 from skimage import img_as_bool
22 from skimage.draw import circle_perimeter 26 from skimage.draw import circle_perimeter
23 from scipy.ndimage.filters import uniform_filter 27 from scipy.ndimage.filters import uniform_filter
24 from scipy.ndimage import gaussian_laplace 28 from scipy.ndimage import gaussian_laplace
25 from os.path import * 29 from os.path import *
523 527
524 return PI2D.getValidOutput() 528 return PI2D.getValidOutput()
525 529
526 530
527 def identifyNumChan(path): 531 def identifyNumChan(path):
528 tiff = tifffile.TiffFile(path) 532
529 shape = tiff.pages[0].shape 533 s = tifffile.TiffFile(path).series[0]
530 numChan=None 534 return s.shape[0] if len(s.shape) > 2 else 1
531 for i, page in enumerate(tiff.pages): 535 # shape = tiff.pages[0].shape
532 if page.shape != shape: 536 # tiff = tifffile.TiffFile(path)
533 numChan = i 537 # for i, page in enumerate(tiff.pages):
534 return numChan 538 # print(page.shape)
535 break 539 # if page.shape != shape:
536 # else: 540 # numChan = i
537 # raise Exception("Did not find any pyramid subresolutions") 541 # return numChan
538 542 # break
539 if not numChan: 543 # else:
540 numChan = len(tiff.pages) 544 # raise Exception("Did not find any pyramid subresolutions")
541 return numChan 545
542 546
543 def getProbMaps(I,dsFactor,modelPath): 547 def getProbMaps(I,dsFactor,modelPath):
544 hsize = int((float(I.shape[0]) * float(0.5))) 548 hsize = int((float(I.shape[0]) * float(0.5)))
545 vsize = int((float(I.shape[1]) * float(0.5))) 549 vsize = int((float(I.shape[1]) * float(0.5)))
546 imagesub = cv2.resize(I,(vsize,hsize),cv2.INTER_NEAREST) 550 imagesub = cv2.resize(I,(vsize,hsize),cv2.INTER_NEAREST)
547 551
548 UNet2D.singleImageInferenceSetup(modelPath, 1) 552 UNet2D.singleImageInferenceSetup(modelPath, 0)
549 553
550 for iSize in range(dsFactor): 554 for iSize in range(dsFactor):
551 hsize = int((float(I.shape[0]) * float(0.5))) 555 hsize = int((float(I.shape[0]) * float(0.5)))
552 vsize = int((float(I.shape[1]) * float(0.5))) 556 vsize = int((float(I.shape[1]) * float(0.5)))
553 I = cv2.resize(I,(vsize,hsize),cv2.INTER_NEAREST) 557 I = cv2.resize(I,(vsize,hsize),cv2.INTER_NEAREST)
555 I = im2double(sk.rescale_intensity(I, in_range=(np.min(I), np.max(I)), out_range=(0, 0.983))) 559 I = im2double(sk.rescale_intensity(I, in_range=(np.min(I), np.max(I)), out_range=(0, 0.983)))
556 probMaps = UNet2D.singleImageInference(I,'accumulate',1) 560 probMaps = UNet2D.singleImageInference(I,'accumulate',1)
557 UNet2D.singleImageInferenceCleanup() 561 UNet2D.singleImageInferenceCleanup()
558 return probMaps 562 return probMaps
559 563
560 def coreSegmenterOutput(I,probMap,initialmask,preBlur,findCenter): 564 def coreSegmenterOutput(I,initialmask,findCenter):
561 hsize = int((float(I.shape[0]) * float(0.1))) 565 hsize = int((float(I.shape[0]) * float(0.1)))
562 vsize = int((float(I.shape[1]) * float(0.1))) 566 vsize = int((float(I.shape[1]) * float(0.1)))
563 nucGF = cv2.resize(I,(vsize,hsize),cv2.INTER_CUBIC) 567 nucGF = cv2.resize(I,(vsize,hsize),cv2.INTER_CUBIC)
564 # Irs = cv2.resize(I,(vsize,hsize),cv2.INTER_CUBIC)
565 # I=I.astype(np.float)
566 # r,c = I.shape
567 # I+=np.random.rand(r,c)*1e-6
568 # c1 = uniform_filter(I, 3, mode='reflect')
569 # c2 = uniform_filter(I*I, 3, mode='reflect')
570 # nucGF = np.sqrt(c2 - c1*c1)*np.sqrt(9./8)
571 # nucGF[np.isnan(nucGF)]=0
572 #active contours 568 #active contours
573 hsize = int(float(nucGF.shape[0])) 569 hsize = int(float(nucGF.shape[0]))
574 vsize = int(float(nucGF.shape[1])) 570 vsize = int(float(nucGF.shape[1]))
575 initialmask = cv2.resize(initialmask,(vsize,hsize),cv2.INTER_NEAREST) 571 initialmask = cv2.resize(initialmask,(vsize,hsize),cv2.INTER_NEAREST)
576 initialmask = dilation(initialmask,disk(15)) >0 572 initialmask = dilation(initialmask,disk(15)) >0
577 573
578 # init=np.argwhere(eroded>0)
579 nucGF = gaussian(nucGF,0.7) 574 nucGF = gaussian(nucGF,0.7)
580 nucGF=nucGF/np.amax(nucGF) 575 nucGF=nucGF/np.amax(nucGF)
581 576
582
583 # initialmask = nucGF>0
584 nuclearMask = morphological_chan_vese(nucGF, 100, init_level_set=initialmask, smoothing=10,lambda1=1.001, lambda2=1) 577 nuclearMask = morphological_chan_vese(nucGF, 100, init_level_set=initialmask, smoothing=10,lambda1=1.001, lambda2=1)
585 578
586 # nuclearMask = chan_vese(nucGF, mu=1.5, lambda1=6, lambda2=1, tol=0.0005, max_iter=2000, dt=15, init_level_set=initialmask, extended_output=True)
587 # nuclearMask = nuclearMask[0]
588
589
590 TMAmask = nuclearMask 579 TMAmask = nuclearMask
591 # nMaskDist =distance_transform_edt(nuclearMask)
592 # fgm = peak_local_max(h_maxima(nMaskDist, 2*preBlur),indices =False)
593 # markers= np.logical_or(erosion(1-nuclearMask,disk(3)),fgm)
594 # TMAmask=watershed(-nMaskDist,label(markers),watershed_line=True)
595 # TMAmask = nuclearMask*(TMAmask>0)
596 TMAmask = remove_small_objects(TMAmask>0,round(TMAmask.shape[0])*round(TMAmask.shape[1])*0.005) 580 TMAmask = remove_small_objects(TMAmask>0,round(TMAmask.shape[0])*round(TMAmask.shape[1])*0.005)
597 TMAlabel = label(TMAmask) 581 TMAlabel = label(TMAmask)
598 # find object closest to center 582 # find object closest to center
599 if findCenter==True: 583 if findCenter==True:
600 584
630 if __name__ == '__main__': 614 if __name__ == '__main__':
631 parser=argparse.ArgumentParser() 615 parser=argparse.ArgumentParser()
632 parser.add_argument("--imagePath") 616 parser.add_argument("--imagePath")
633 parser.add_argument("--outputPath") 617 parser.add_argument("--outputPath")
634 parser.add_argument("--maskPath") 618 parser.add_argument("--maskPath")
635 parser.add_argument("--downsampleFactor",type = int, default = 5) 619 parser.add_argument("--tissue", action='store_true')
620 parser.add_argument("--downsampleFactor", type=int, default = 5)
636 parser.add_argument("--channel",type = int, default = 0) 621 parser.add_argument("--channel",type = int, default = 0)
637 parser.add_argument("--buffer",type = float, default = 2) 622 parser.add_argument("--buffer",type = float, default = 2)
638 parser.add_argument("--outputChan", type=int, nargs = '+', default=[-1]) 623 parser.add_argument("--outputChan", type=int, nargs = '+', default=[-1])
639 parser.add_argument("--sensitivity",type = float, default=0.3) 624 parser.add_argument("--sensitivity",type = float, default=0.3)
640 parser.add_argument("--useGrid",action='store_true') 625 parser.add_argument("--useGrid",action='store_true')
641 parser.add_argument("--cluster",action='store_true') 626 parser.add_argument("--cluster",action='store_true')
642 args = parser.parse_args() 627 args = parser.parse_args()
643 628
644 outputPath = args.outputPath 629 outputPath = args.outputPath
645 imagePath = args.imagePath 630 imagePath = args.imagePath
646 sensitivity = args.sensitivity 631 sensitivity = args.sensitivity
647 #scriptPath = os.path.dirname(os.path.realpath(__file__))
648 #modelPath = os.path.join(scriptPath, 'TFModel - 3class 16 kernels 5ks 2 layers')
649 #modelPath = 'D:\\LSP\\Coreograph\\model-4layersMaskAug20'
650 scriptPath = os.path.dirname(os.path.realpath(__file__)) 632 scriptPath = os.path.dirname(os.path.realpath(__file__))
651 modelPath = os.path.join(scriptPath, 'model') 633 modelPath = os.path.join(scriptPath, 'model')
652 # outputPath = 'D:\\LSP\\cycif\\testsets\\exemplar-002\\dearrayPython' ############
653 maskOutputPath = os.path.join(outputPath, 'masks') 634 maskOutputPath = os.path.join(outputPath, 'masks')
654 # imagePath = 'D:\\LSP\\cycif\\testsets\\exemplar-002\\registration\\exemplar-002.ome.tif'###########
655 # imagePath = 'Y:\\sorger\\data\\RareCyte\\Connor\\TMAs\\CAJ_TMA11_13\\original_data\\TMA11\\registration\\TMA11.ome.tif'
656 # imagePath = 'Y:\\sorger\\data\\RareCyte\\Connor\\TMAs\\Z124_TMA20_22\\TMA22\\registration\\TMA22.ome.tif'
657 # classProbsPath = 'D:\\unetcoreograph.tif'
658 # imagePath = 'Y:\\sorger\\data\\RareCyte\\Connor\\Z155_PTCL\\TMA_552\\registration\\TMA_552.ome.tif'
659 # classProbsPath = 'Y:\\sorger\\data\\RareCyte\\Connor\\Z155_PTCL\\TMA_552\\probMapCore\\TMA_552_CorePM_1.tif'
660 # imagePath = 'Y:\\sorger\\data\\RareCyte\\Zoltan\\Z112_TMA17_19\\190403_ashlar\\TMA17_1092.ome.tif'
661 # classProbsPath = 'Z:\\IDAC\\Clarence\\LSP\\CyCIF\\TMA\\probMapCore\\1new_CorePM_1.tif'
662 # imagePath = 'Y:\\sorger\\data\\RareCyte\\ANNIINA\\Julia\\2018\\TMA6\\julia_tma6.ome.tif'
663 # classProbsPath = 'Z:\\IDAC\\Clarence\\LSP\\CyCIF\\TMA\\probMapCore\\3new_CorePM_1.tif'
664 635
665 636
666 # if not os.path.exists(outputPath): 637 # if not os.path.exists(outputPath):
667 # os.makedirs(outputPath) 638 # os.makedirs(outputPath)
668 # else: 639 # else:
669 # shutil.rmtree(outputPath) 640 # shutil.rmtree(outputPath)
670 if not os.path.exists(maskOutputPath): 641 if not os.path.exists(maskOutputPath):
671 os.makedirs(maskOutputPath) 642 os.makedirs(maskOutputPath)
672 643 print(
673 644 'WARNING! IF USING FOR TISSUE SPLITTING, IT IS ADVISED TO SET --downsampleFactor TO HIGHER THAN DEFAULT OF 5')
674 channel = args.channel 645 channel = args.channel
675 dsFactor = 1/(2**args.downsampleFactor) 646 dsFactor = 1/(2**args.downsampleFactor)
676 # I = tifffile.imread(imagePath, key=channel)
677 I = skio.imread(imagePath, img_num=channel) 647 I = skio.imread(imagePath, img_num=channel)
678
679 imagesub = resize(I,(int((float(I.shape[0]) * dsFactor)),int((float(I.shape[1]) * dsFactor)))) 648 imagesub = resize(I,(int((float(I.shape[0]) * dsFactor)),int((float(I.shape[1]) * dsFactor))))
680 numChan = identifyNumChan(imagePath) 649 numChan = identifyNumChan(imagePath)
681 650
682 outputChan = args.outputChan 651 outputChan = args.outputChan
683 if len(outputChan)==1: 652 if len(outputChan)==1:
684 if outputChan[0]==-1: 653 if outputChan[0]==-1:
685 outputChan = [0, numChan-1] 654 outputChan = [0, numChan-1]
686 else: 655 else:
687 outputChan.append(outputChan[0]) 656 outputChan.append(outputChan[0])
688 657 classProbs = getProbMaps(I, args.downsampleFactor, modelPath)
689 classProbs = getProbMaps(I,args.downsampleFactor,modelPath) 658
690 # classProbs = tifffile.imread(classProbsPath,key=0) 659 if not args.tissue:
691 preMask = gaussian(np.uint8(classProbs*255),1)>0.8 660 print('TMA mode selected')
692 661 preMask = gaussian(np.uint8(classProbs*255),1)>0.8
693 P = regionprops(label(preMask),cache=False) 662
694 area = [ele.area for ele in P] 663 P = regionprops(label(preMask),cache=False)
695 print(str(len(P)) + ' cores detected!') 664 area = [ele.area for ele in P]
696 if len(P) <3: 665 if len(P) <3:
697 medArea = np.median(area) 666 medArea = np.median(area)
698 maxArea = np.percentile(area,99) 667 maxArea = np.percentile(area,99)
668 else:
669 count=0
670 labelpreMask = np.zeros(preMask.shape,dtype=np.uint32)
671 for props in P:
672 count += 1
673 yi = props.coords[:, 0]
674 xi = props.coords[:, 1]
675 labelpreMask[yi, xi] = count
676 P=regionprops(labelpreMask)
677 area = [ele.area for ele in P]
678 medArea = np.median(area)
679 maxArea = np.percentile(area,99)
680 preMask = remove_small_objects(preMask,0.2*medArea)
681 coreRad = round(np.sqrt(medArea/np.pi))
682 estCoreDiam = round(np.sqrt(maxArea/np.pi)*1.2*args.buffer)
683
684 #preprocessing
685 fgFiltered = blob_log(preMask,coreRad*0.6,threshold=sensitivity)
686 Imax = np.zeros(preMask.shape,dtype=np.uint8)
687 for iSpot in range(fgFiltered.shape[0]):
688 yi = np.uint32(round(fgFiltered[iSpot, 0]))
689 xi = np.uint32(round(fgFiltered[iSpot, 1]))
690 Imax[yi, xi] = 1
691 Imax = Imax*preMask
692 Idist = distance_transform_edt(1-Imax)
693 markers = label(Imax)
694 coreLabel = watershed(Idist,markers,watershed_line=True,mask = preMask)
695 P = regionprops(coreLabel)
696 centroids = np.array([ele.centroid for ele in P]) / dsFactor
697 np.savetxt(outputPath + os.path.sep + 'centroidsY-X.txt', np.asarray(centroids), fmt='%10.5f')
698 numCores = len(centroids)
699 print(str(numCores) + ' cores detected!')
700 estCoreDiamX = np.ones(numCores) * estCoreDiam / dsFactor
701 estCoreDiamY = np.ones(numCores) * estCoreDiam / dsFactor
699 else: 702 else:
700 count=0 703 print('Tissue mode selected')
701 labelpreMask = np.zeros(preMask.shape,dtype=np.uint32) 704 imageblur = 5
702 for props in P: 705 Iblur = gaussian(np.uint8(255*classProbs), imageblur)
703 count += 1 706 coreMask = binary_fill_holes(binary_closing(Iblur > threshold_otsu(Iblur), np.ones((imageblur*2,imageblur*2))))
704 yi = props.coords[:, 0] 707 coreMask = remove_small_objects(coreMask, min_size=0.001 * coreMask.shape[0] * coreMask.shape[1])
705 xi = props.coords[:, 1] 708
706 labelpreMask[yi, xi] = count 709 ## watershed
707 P=regionprops(labelpreMask) 710 Idist = distance_transform_edt(coreMask)
708 area = [ele.area for ele in P] 711 markers = peak_local_max(h_maxima(Idist,20),indices=False)
709 medArea = np.median(area) 712 markers = label(markers).astype(np.int8)
710 maxArea = np.percentile(area,99) 713 coreLabel = watershed(-Idist, markers, watershed_line=True,mask = coreMask)
711 preMask = remove_small_objects(preMask,0.2*medArea) 714
712 coreRad = round(np.sqrt(medArea/np.pi)) 715 P = regionprops(coreLabel)
713 estCoreDiam = round(np.sqrt(maxArea/np.pi)*1.2*args.buffer) 716 centroids = np.array([ele.centroid for ele in P]) / dsFactor
714 717 np.savetxt(outputPath + os.path.sep + 'centroidsY-X.txt', np.asarray(centroids), fmt='%10.5f')
715 #preprocessing 718 numCores = len(centroids)
716 fgFiltered = blob_log(preMask,coreRad*0.6,threshold=sensitivity) 719 print(str(numCores) + ' tissues detected!')
717 Imax = np.zeros(preMask.shape,dtype=np.uint8) 720 estCoreDiamX = np.array([(ele.bbox[3]-ele.bbox[1])*1.1 for ele in P]) / dsFactor
718 for iSpot in range(fgFiltered.shape[0]): 721 estCoreDiamY = np.array([(ele.bbox[2]-ele.bbox[0])*1.1 for ele in P]) / dsFactor
719 yi = np.uint32(round(fgFiltered[iSpot, 0]))
720 xi = np.uint32(round(fgFiltered[iSpot, 1]))
721 Imax[yi, xi] = 1
722 Imax = Imax*preMask
723 Idist = distance_transform_edt(1-Imax)
724 markers = label(Imax)
725 coreLabel = watershed(Idist,markers,watershed_line=True,mask = preMask)
726 P = regionprops(coreLabel)
727 centroids = np.array([ele.centroid for ele in P])/dsFactor
728 numCores = len(centroids)
729 estCoreDiamX = np.ones(numCores)*estCoreDiam/dsFactor
730 estCoreDiamY = np.ones(numCores)*estCoreDiam/dsFactor
731 722
732 if numCores ==0 & args.cluster: 723 if numCores ==0 & args.cluster:
733 print('No cores detected. Try adjusting the downsample factor') 724 print('No cores detected. Try adjusting the downsample factor')
734 sys.exit(255) 725 sys.exit(255)
735 726
736 singleMaskTMA = np.zeros(imagesub.shape) 727 singleMaskTMA = np.zeros(imagesub.shape)
737 maskTMA = np.zeros(imagesub.shape) 728 maskTMA = np.zeros(imagesub.shape)
738 bbox = [None] * numCores 729 bbox = [None] * numCores
739 730 imagesub = imagesub/np.percentile(imagesub,99.9)
740 731 imagesub = (imagesub * 255).round().astype(np.uint8)
732 imagesub = gray2rgb(imagesub)
741 x=np.zeros(numCores) 733 x=np.zeros(numCores)
742 xLim=np.zeros(numCores) 734 xLim=np.zeros(numCores)
743 y=np.zeros(numCores) 735 y=np.zeros(numCores)
744 yLim=np.zeros(numCores) 736 yLim=np.zeros(numCores)
745 737
759 yLim[iCore] = I.shape[0] 751 yLim[iCore] = I.shape[0]
760 if y[iCore]<1: 752 if y[iCore]<1:
761 y[iCore]=1 753 y[iCore]=1
762 754
763 bbox[iCore] = [round(x[iCore]), round(y[iCore]), round(xLim[iCore]), round(yLim[iCore])] 755 bbox[iCore] = [round(x[iCore]), round(y[iCore]), round(xLim[iCore]), round(yLim[iCore])]
764 756 coreStack = np.zeros((outputChan[1]-outputChan[0]+1,np.int(round(yLim[iCore])-round(y[iCore])-1),np.int(round(xLim[iCore])-round(x[iCore])-1)),dtype='uint16')
757
765 for iChan in range(outputChan[0],outputChan[1]+1): 758 for iChan in range(outputChan[0],outputChan[1]+1):
766 with pytiff.Tiff(imagePath, "r", encoding='utf-8') as handle: 759 with pytiff.Tiff(imagePath, "r", encoding='utf-8') as handle:
767 handle.set_page(iChan) 760 handle.set_page(iChan)
768 coreStack= handle[np.uint32(bbox[iCore][1]):np.uint32(bbox[iCore][3]-1), np.uint32(bbox[iCore][0]):np.uint32(bbox[iCore][2]-1)] 761 coreStack[iChan,:,:] =handle[np.uint32(bbox[iCore][1]):np.uint32(bbox[iCore][3]-1), np.uint32(bbox[iCore][0]):np.uint32(bbox[iCore][2]-1)]
769 skio.imsave(outputPath + os.path.sep + str(iCore+1) + '.tif',coreStack,append=True) 762
770 763 skio.imsave(outputPath + os.path.sep + str(iCore+1) + '.tif',np.uint16(coreStack),imagej=True,bigtiff=True)
771 with pytiff.Tiff(imagePath, "r", encoding='utf-8') as handle: 764 with pytiff.Tiff(imagePath, "r", encoding='utf-8') as handle:
772 handle.set_page(args.channel) 765 handle.set_page(args.channel)
773 coreSlice= handle[np.uint32(bbox[iCore][1]):np.uint32(bbox[iCore][3]-1), np.uint32(bbox[iCore][0]):np.uint32(bbox[iCore][2]-1)] 766 coreSlice= handle[np.uint32(bbox[iCore][1]):np.uint32(bbox[iCore][3]-1), np.uint32(bbox[iCore][0]):np.uint32(bbox[iCore][2]-1)]
774 767
775 core = (coreLabel ==(iCore+1)) 768 core = (coreLabel ==(iCore+1))
776 initialmask = core[np.uint32(y[iCore]*dsFactor):np.uint32(yLim[iCore]*dsFactor),np.uint32(x[iCore]*dsFactor):np.uint32(xLim[iCore]*dsFactor)] 769 initialmask = core[np.uint32(y[iCore] * dsFactor):np.uint32(yLim[iCore] * dsFactor),
777 initialmask = resize(initialmask,size(coreSlice),cv2.INTER_NEAREST) 770 np.uint32(x[iCore] * dsFactor):np.uint32(xLim[iCore] * dsFactor)]
778 771 if not args.tissue:
779 singleProbMap = classProbs[np.uint32(y[iCore]*dsFactor):np.uint32(yLim[iCore]*dsFactor),np.uint32(x[iCore]*dsFactor):np.uint32(xLim[iCore]*dsFactor)] 772 initialmask = resize(initialmask,size(coreSlice),cv2.INTER_NEAREST)
780 singleProbMap = resize(np.uint8(255*singleProbMap),size(coreSlice),cv2.INTER_NEAREST) 773
781 TMAmask = coreSegmenterOutput(coreSlice,singleProbMap,initialmask,coreRad/20,False) 774 singleProbMap = classProbs[np.uint32(y[iCore]*dsFactor):np.uint32(yLim[iCore]*dsFactor),np.uint32(x[iCore]*dsFactor):np.uint32(xLim[iCore]*dsFactor)]
775 singleProbMap = resize(np.uint8(255*singleProbMap),size(coreSlice),cv2.INTER_NEAREST)
776 TMAmask = coreSegmenterOutput(coreSlice,initialmask,False)
777 else:
778 Irs = resize(coreSlice,(int((float(coreSlice.shape[0]) * 0.25)),int((float(coreSlice.shape[1]) * 0.25))))
779 TMAmask = coreSegmenterOutput(Irs, np.uint8(initialmask), False)
780
782 if np.sum(TMAmask)==0: 781 if np.sum(TMAmask)==0:
783 TMAmask = np.ones(TMAmask.shape) 782 TMAmask = np.ones(TMAmask.shape)
784 vsize = int(float(coreSlice.shape[0])) 783 vsize = int(float(coreSlice.shape[0]))
785 hsize = int(float(coreSlice.shape[1])) 784 hsize = int(float(coreSlice.shape[1]))
786 masksub = resize(resize(TMAmask,(vsize,hsize),cv2.INTER_NEAREST),(int((float(coreSlice.shape[0])*dsFactor)),int((float(coreSlice.shape[1])*dsFactor))),cv2.INTER_NEAREST) 785 masksub = resize(resize(TMAmask,(vsize,hsize),cv2.INTER_NEAREST),(int((float(coreSlice.shape[0])*dsFactor)),int((float(coreSlice.shape[1])*dsFactor))),cv2.INTER_NEAREST)
787 singleMaskTMA[int(y[iCore]*dsFactor):int(y[iCore]*dsFactor)+masksub.shape[0],int(x[iCore]*dsFactor):int(x[iCore]*dsFactor)+masksub.shape[1]]=masksub 786 singleMaskTMA[int(y[iCore]*dsFactor):int(y[iCore]*dsFactor)+masksub.shape[0],int(x[iCore]*dsFactor):int(x[iCore]*dsFactor)+masksub.shape[1]]=masksub
788 maskTMA = maskTMA + resize(singleMaskTMA,maskTMA.shape,cv2.INTER_NEAREST) 787 maskTMA = maskTMA + resize(singleMaskTMA,maskTMA.shape,cv2.INTER_NEAREST)
789 cv2.putText(imagesub, str(iCore+1), (int(P[iCore].centroid[1]),int(P[iCore].centroid[0])), 0, 0.5, (np.amax(imagesub), np.amax(imagesub), np.amax(imagesub)), 1, cv2.LINE_AA) 788
789 cv2.putText(imagesub, str(iCore+1), (int(P[iCore].centroid[1]),int(P[iCore].centroid[0])), 0, 0.5, (0,255,0), 1, cv2.LINE_AA)
790 790
791 skio.imsave(maskOutputPath + os.path.sep + str(iCore+1) + '_mask.tif',np.uint8(TMAmask)) 791 skio.imsave(maskOutputPath + os.path.sep + str(iCore+1) + '_mask.tif',np.uint8(TMAmask))
792 print('Segmented core ' + str(iCore+1)) 792 print('Segmented core/tissue ' + str(iCore+1))
793 793
794 boundaries = find_boundaries(maskTMA) 794 boundaries = find_boundaries(maskTMA)
795 imagesub = imagesub/np.percentile(imagesub,99.9) 795 imagesub[boundaries==1] = 255
796 imagesub[boundaries==1] = 1 796 skio.imsave(outputPath + os.path.sep + 'TMA_MAP.tif' ,imagesub)
797 skio.imsave(outputPath + os.path.sep + 'TMA_MAP.tif' ,np.uint8(imagesub*255)) 797 print('Segmented all cores/tissues!')
798 print('Segmented all cores!')
799
800 798
801 #restore GPU to 0 799 #restore GPU to 0
802 #image load using tifffile 800 #image load using tifffile