diff --git a/analysis/database.py b/analysis/database.py index 1e909a74577ce423fb2fa2feef958156dc6a5a8d..b37bb1dceafbdf8ed14ae2fe3ebbd634dc7e2c49 100644 --- a/analysis/database.py +++ b/analysis/database.py @@ -44,7 +44,7 @@ class DataBaseWindow(QtWidgets.QMainWindow): self.path = os.path.join(Path.home(), 'gepard', 'databases') self.importPath = self.path if not os.path.exists(self.path): - os.mkdir(self.path) + os.makedirs(self.path) self.activeDatabase = None self.activeSpectrum = None self.activeSpectrumName = None diff --git a/dataset.py b/dataset.py index f11f44762f7daa1b1bd614acd5c181738f49fa4e..37728185aaac02fa12866e1d69232bf4d0ea649b 100644 --- a/dataset.py +++ b/dataset.py @@ -1,387 +1,412 @@ -# -*- coding: utf-8 -*- -""" -GEPARD - Gepard-Enabled PARticle Detection -Copyright (C) 2018 Lars Bittrich and Josef Brandt, Leibniz-Institut für -Polymerforschung Dresden e. V. - -This program is free software: you can redistribute it and/or modify -it under the terms of the GNU General Public License as published by -the Free Software Foundation, either version 3 of the License, or -(at your option) any later version. - -This program is distributed in the hope that it will be useful, -but WITHOUT ANY WARRANTY; without even the implied warranty of -MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -GNU General Public License for more details. - -You should have received a copy of the GNU General Public License -along with this program, see COPYING. -If not, see . -""" -import os -import pickle -import numpy as np -import cv2 +# -*- coding: utf-8 -*- +""" +GEPARD - Gepard-Enabled PARticle Detection +Copyright (C) 2018 Lars Bittrich and Josef Brandt, Leibniz-Institut für +Polymerforschung Dresden e. V. + +This program is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 3 of the License, or +(at your option) any later version. + +This program is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with this program, see COPYING. +If not, see . +""" +import os +import pickle +import numpy as np +import cv2 from .helperfunctions import cv2imread_fix, cv2imwrite_fix -from copy import copy - -currentversion = 2 - -def loadData(fname): - retds = None - with open(fname, "rb") as fp: - ds = pickle.load(fp) - ds.fname = fname - ds.readin = True - ds.updatePath() - retds = DataSet(fname) - retds.version = 0 - retds.__dict__.update(ds.__dict__) - if retds.version < currentversion: - retds.legacyConversion() - elif retds.zvalimg=="saved": - retds.loadZvalImg() - return retds - -def saveData(dataset, fname): - with open(fname, "wb") as fp: - # zvalimg is rather large and thus it is saved separately in a tif file - # only onces after its creation - zvalimg = dataset.zvalimg - if zvalimg is not None: - dataset.zvalimg = "saved" - pickle.dump(dataset, fp, protocol=-1) - dataset.zvalimg = zvalimg - -def arrayCompare(a1, a2): - if a1.shape!=a2.shape: - return False - if a1.dtype!=np.float32 and a1.dtype!=np.float64: - return np.all(a1==a2) - ind = np.isnan(a1) - if not np.any(ind): - return np.all(a1==a2) - return np.all(a1[~ind]==a2[~ind]) - -def listCompare(l1, l2): - if len(l1)!=len(l2): - return False - for l1i, l2i in zip(l1, l2): - if isinstance(l1i, np.ndarray): - if not isinstance(l2i, np.ndarray) or not arrayCompare(l1i, l2i): - return False - elif isinstance(l1i, (list, tuple)): - if not isinstance(l2i, (list, tuple)) or not listCompare(l1i, l2i): - return False - elif l1i!=l2i and ((~np.isnan(l1i)) or (~np.isnan(l2i))): - return False - return True - -def recursiveDictCompare(d1, d2): - for key in d1: - if not key in d2: - print("key missing in d2:", key, flush=True) - return False - a = d1[key] - b = d2[key] - print(key, type(a), type(b), flush=True) - if isinstance(a, np.ndarray): - if not isinstance(b, np.ndarray) or not arrayCompare(a, b): - print("data is different!", a, b) - return False - elif isinstance(a, dict): - if not isinstance(b, dict): - print("data is different!", a, b) - return False - if not recursiveDictCompare(a, b): - return False - elif isinstance(a, (list, tuple)): - if not isinstance(b, (list, tuple)) or not listCompare(a, b): - print("data is different!", a, b) - return False - elif a != b: - if (a is not None) and (b is not None): - print("data is different!", a, b) - return False - return True - -class DataSet(object): - def __init__(self, fname, newProject=False): - self.fname = fname - # parameters specifically for optical scan - self.version = currentversion - self.lastpos = None - self.maxdim = None - self.pixelscale_df = None # µm / pixel --> scale of DARK FIELD camera (used for image stitching) - self.pixelscale_bf = None # µm / pixel of DARK FIELD camera (set to same as bright field, if both use the same camera) - self.imagedim_bf = None # width, height, angle of BRIGHT FIELD camera - self.imagedim_df = None # width, height, angle of DARK FIELD camera (set to same as bright field, if both use the same camera) - self.imagescanMode = 'df' #was the fullimage acquired in dark- or brightfield? - self.fitpoints = [] # manually adjusted positions aquired to define the specimen geometry - self.fitindices = [] # which of the five positions in the ui are already known - self.boundary = [] # scan boundary computed by a circle around the fitpoints + manual adjustments - self.grid = [] # scan grid positions for optical scan - self.zpositions = [] # z-positions for optical scan - self.heightmap = None - self.zvalimg = None - - # parameters specifically for raman scan - self.pshift = None # shift of raman scan position relative to image center - self.coordOffset = [0, 0] #offset of entire coordinate system - self.seedpoints = np.array([]) - self.seeddeletepoints = np.array([]) - self.detectParams = {'points': np.array([[50,0],[100,200],[200,255]]), - 'contrastcurve': True, - 'blurRadius': 9, - 'threshold': 0.2, - 'maxholebrightness': 0.5, - 'erodeconvexdefects': 0, - 'minparticlearea': 20, - 'minparticledistance': 20, - 'measurefrac': 1, - 'compactness': 0.1, - 'seedRad': 3} - - self.ramanpoints = [] - self.particlecontours = [] - self.particlestats = [] - self.ramanscansortindex = None - self.ramanscandone = False - - self.results = {'polymers': None, - 'hqis': None, - 'additives': None, - 'additive_hqis': None} - - self.resultParams = {'minHQI': None, - 'compHQI': None} - self.spectraPath = None - self.particles2spectra = None #links idParticle to corresponding idSpectra (i.e., first measured particle (ID=0) is linked to spectra indices 0 and 1) - self.colorSeed = 'default' - self.resultsUploadedToSQL = [] - - self.readin = True # a value that is always set to True at loadData - # and mark that the coordinate system might be changed in the meantime - self.mode = "prepare" - if newProject: - self.fname = self.newProject(fname) - self.updatePath() - - def __eq__(self, other): - return recursiveDictCompare(self.__dict__, other.__dict__) - - def getPixelScale(self, mode=None): - if mode is None: - mode = self.imagescanMode - return (self.pixelscale_df if mode == 'df' else self.pixelscale_bf) - - def saveZvalImg(self): - if self.zvalimg is not None: - cv2imwrite_fix(self.getZvalImageName(), self.zvalimg) - - def loadZvalImg(self): - if os.path.exists(self.getZvalImageName()): - self.zvalimg = cv2imread_fix(self.getZvalImageName(), cv2.IMREAD_GRAYSCALE) - - def legacyConversion(self, recreatefullimage=False): - if self.version==0: - print("Converting legacy version 0 to 1") - print("This may take some time") - - # local imports as these functions are only needed for the rare occasion of legacy conversion - from opticalscan import loadAndPasteImage - - # try to load png and check for detection contours - recreatefullimage = recreatefullimage or not os.path.exists(self.getLegacyImageName()) - if not recreatefullimage: - img = cv2imread_fix(self.getLegacyImageName()) - Nc = len(self.particlecontours) - if Nc>0: - contour = self.particlecontours[Nc//2] - contpixels = img[contour[:,0,1],contour[:,0,0]] - if np.all(contpixels[:,1]==255) and np.all(contpixels[:,2]==0) \ - and np.all(contpixels[:,0]==0): - recreatefullimage = True - if not recreatefullimage: - cv2imwrite_fix(self.getImageName(), img) - del img - - if recreatefullimage: - print("recreating fullimage from grid data") - imgdata = None - zvalimg = None - Ngrid = len(self.grid) - - width, height, rotationvalue = self.imagedim_df - p0, p1 = self.maxdim[:2], self.maxdim[2:] - for i in range(Ngrid): - print(f"Processing image {i+1} of {Ngrid}") - names = [] - for k in range(len(self.zpositions)): - names.append(os.path.join(self.getScanPath(), f"image_{i}_{k}.bmp")) - p = self.grid[i] - imgdata, zvalimg = loadAndPasteImage(names, imgdata, zvalimg, width, - height, rotationvalue, p0, p1, p) - self.zvalimg = zvalimg - cv2imwrite_fix(self.getImageName(), cv2.cvtColor(imgdata, cv2.COLOR_RGB2BGR)) - del imgdata - self.saveZvalImg() - if "particleimgs" in self.__dict__: - del self.particleimgs - - self.version = 1 - - - if self.version == 1: - print("Converting legacy version 1 to 2") - if hasattr(self, 'pixelscale'): - print('pixelscale was', self.pixelscale) - self.pixelscale_bf = self.pixelscale - self.pixelscale_df = self.pixelscale - del self.pixelscale - - if hasattr(self, 'imagedim'): - self.imagedim_bf = self.imagedim - self.imagedim_df = self.imagedim - del self.imagedim - - self.version = 2 - - # add later conversion for higher version numbers here - - def getSubImage(self, img, index, draw=True): - contour = self.particlecontours[index] - x0, x1 = contour[:,0,0].min(), contour[:,0,0].max() - y0, y1 = contour[:,0,1].min(), contour[:,0,1].max() - subimg = img[y0:y1+1,x0:x1+1].copy() - if draw: - cv2.drawContours(subimg, [contour], -1, (0,255,0), 1) - return subimg - - def getZval(self, pixelpos): - assert self.zvalimg is not None - zp = self.zvalimg[round(pixelpos[1]), round(pixelpos[0])] - z0, z1 = self.zpositions.min(), self.zpositions.max() - return zp/255.*(z1-z0) + z0 - - def mapHeight(self, x, y): - assert not self.readin - assert self.heightmap is not None - return self.heightmap[0]*x + self.heightmap[1]*y + self.heightmap[2] - - def mapToPixel(self, p, mode='df', force=False): - if not force: - assert not self.readin - p0 = copy(self.lastpos) - - if mode == 'df': - p0[0] -= self.imagedim_df[0]/2 - p0[1] += self.imagedim_df[1]/2 - return (p[0] - p0[0])/self.pixelscale_df, (p0[1] - p[1])/self.pixelscale_df - - elif mode == 'bf': - p0[0] -= self.imagedim_bf[0]/2 - p0[1] += self.imagedim_bf[1]/2 - return (p[0] - p0[0])/self.pixelscale_bf, (p0[1] - p[1])/self.pixelscale_bf - else: - print('mapToPixelMode not understood') - return - - def mapToLength(self, pixelpos, mode='df', force=False): - if not force: - assert not self.readin - p0 = copy(self.lastpos) - p0[0] += self.coordOffset[0] - p0[1] += self.coordOffset[1] - - if mode == 'df': - p0[0] -= self.imagedim_df[0]/2 - p0[1] += self.imagedim_df[1]/2 - return (pixelpos[0]*self.pixelscale_df + p0[0]), (p0[1] - pixelpos[1]*self.pixelscale_df) - elif mode == 'bf': - p0[0] -= self.imagedim_bf[0]/2 - p0[1] += self.imagedim_bf[1]/2 - return (pixelpos[0]*self.pixelscale_bf + p0[0]), (p0[1] - pixelpos[1]*self.pixelscale_bf) - else: - raise ValueError(f'mapToLength mode: {mode} not understood') - - def mapToLengthRaman(self, pixelpos, microscopeMode='df', noz=False): - p0x, p0y = self.mapToLength(pixelpos, mode = microscopeMode) - x, y = p0x + self.pshift[0], p0y + self.pshift[1] - z = None - if not noz: - z = self.mapHeight(x, y) - z += self.getZval(pixelpos) - return x, y, z - - def newProject(self, fname): - path = os.path.split(fname)[0] - name = os.path.splitext(os.path.basename(fname))[0] - newpath = os.path.join(path, name) - fname = os.path.join(newpath, name + ".pkl") - if not os.path.exists(newpath): - os.mkdir(newpath) # for new projects a directory will be created - elif os.path.exists(fname): # if this project is already there, load it instead - self.__dict__.update(loadData(fname).__dict__) - return fname - - def getScanPath(self): - scandir = os.path.join(self.path, "scanimages") - if not os.path.exists(scandir): - os.mkdir(scandir) - return scandir - - def updatePath(self): - self.path = os.path.split(self.fname)[0] - self.name = os.path.splitext(os.path.basename(self.fname))[0] - - def getImageName(self): - return os.path.join(self.path, 'fullimage.tif') - - def getZvalImageName(self): - return os.path.join(self.path, "zvalues.tif") - - def getLegacyImageName(self): - return os.path.join(self.path, "fullimage.png") - - def getLegacyDetectImageName(self): - return os.path.join(self.path, "detectimage.png") - - def getDetectImageName(self): - raise NotImplementedError("No longer implemented due to change in API") - - def getTmpImageName(self): - return os.path.join(self.path, "tmp.bmp") - - def saveParticleData(self): - print('Not saving ParticleData into text file...:\nThe current output format might be wrong, if multiple spectra per particle are present...') -# if len(self.ramanscansortindex)>0: -# data = [] -# pixelscale = (self.pixelscale_df if self.imagescanMode == 'df' else self.pixelscale_bf) -# for i in self.ramanscansortindex: -# data.append(list(self.ramanpoints[i])+list(self.particlestats[i])) -# data = np.array(data) -# data[:,0], data[:,1], z = self.mapToLengthRaman((data[:,0], data[:,1]), microscopeMode=self.imagescanMode, noz=True) -# data[:,2:7] *= pixelscale -# header = "x [µm], y [µm], length [µm], height [µm], length_ellipse [µm], height_ellipse [µm]" -# if data.shape[1]>6: -# header = header + ", area [µm^2]" -# data[:,6] *= pixelscale -# np.savetxt(os.path.join(self.path, "particledata.txt"), data, -# header=header) - - def save(self): - saveData(self, self.fname) - - def saveBackup(self): - inc = 0 - while True: - directory = os.path.dirname(self.fname) - filename = self.name + '_backup_' + str(inc) + '.pkl' - path = os.path.join(directory, filename) - if os.path.exists(path): - inc += 1 - else: - saveData(self, path) - return filename - +from copy import copy + +currentversion = 2 + +def loadData(fname): + retds = None + with open(fname, "rb") as fp: + ds = pickle.load(fp) + ds.fname = fname + ds.readin = True + ds.updatePath() + retds = DataSet(fname) + retds.version = 0 + retds.__dict__.update(ds.__dict__) + if retds.version < currentversion: + retds.legacyConversion() + elif retds.zvalimg=="saved": + retds.loadZvalImg() + return retds + +def saveData(dataset, fname): + with open(fname, "wb") as fp: + # zvalimg is rather large and thus it is saved separately in a tif file + # only onces after its creation + zvalimg = dataset.zvalimg + if zvalimg is not None: + dataset.zvalimg = "saved" + pickle.dump(dataset, fp, protocol=-1) + dataset.zvalimg = zvalimg + +def arrayCompare(a1, a2): + if a1.shape!=a2.shape: + return False + if a1.dtype!=np.float32 and a1.dtype!=np.float64: + return np.all(a1==a2) + ind = np.isnan(a1) + if not np.any(ind): + return np.all(a1==a2) + return np.all(a1[~ind]==a2[~ind]) + +def listCompare(l1, l2): + if len(l1)!=len(l2): + return False + for l1i, l2i in zip(l1, l2): + if isinstance(l1i, np.ndarray): + if not isinstance(l2i, np.ndarray) or not arrayCompare(l1i, l2i): + return False + elif isinstance(l1i, (list, tuple)): + if not isinstance(l2i, (list, tuple)) or not listCompare(l1i, l2i): + return False + elif l1i!=l2i and ((~np.isnan(l1i)) or (~np.isnan(l2i))): + return False + return True + +def recursiveDictCompare(d1, d2): + for key in d1: + if not key in d2: + print("key missing in d2:", key, flush=True) + return False + a = d1[key] + b = d2[key] + print(key, type(a), type(b), flush=True) + if isinstance(a, np.ndarray): + if not isinstance(b, np.ndarray) or not arrayCompare(a, b): + print("data is different!", a, b) + return False + elif isinstance(a, dict): + if not isinstance(b, dict): + print("data is different!", a, b) + return False + if not recursiveDictCompare(a, b): + return False + elif isinstance(a, (list, tuple)): + if not isinstance(b, (list, tuple)) or not listCompare(a, b): + print("data is different!", a, b) + return False + elif a != b: + if (a is not None) and (b is not None): + print("data is different!", a, b) + return False + return True + +class DataSet(object): + def __init__(self, fname, newProject=False): + self.fname = fname + # parameters specifically for optical scan + self.version = currentversion + self.lastpos = None + self.maxdim = None + self.pixelscale_df = None # µm / pixel --> scale of DARK FIELD camera (used for image stitching) + self.pixelscale_bf = None # µm / pixel of DARK FIELD camera (set to same as bright field, if both use the same camera) + self.imagedim_bf = None # width, height, angle of BRIGHT FIELD camera + self.imagedim_df = None # width, height, angle of DARK FIELD camera (set to same as bright field, if both use the same camera) + self.imagescanMode = 'df' #was the fullimage acquired in dark- or brightfield? + self.fitpoints = [] # manually adjusted positions aquired to define the specimen geometry + self.fitindices = [] # which of the five positions in the ui are already known + self.boundary = [] # scan boundary computed by a circle around the fitpoints + manual adjustments + self.grid = [] # scan grid positions for optical scan + self.zpositions = [] # z-positions for optical scan + self.heightmap = None + self.zvalimg = None + self.coordinatetransform = None # if imported form extern source coordinate system may be rotated + self.signx = 1. + self.signy = -1. + + # parameters specifically for raman scan + self.pshift = None # shift of raman scan position relative to image center + self.coordOffset = [0, 0] #offset of entire coordinate system + self.seedpoints = np.array([]) + self.seeddeletepoints = np.array([]) + self.detectParams = {'points': np.array([[50,0],[100,200],[200,255]]), + 'contrastcurve': True, + 'blurRadius': 9, + 'threshold': 0.2, + 'maxholebrightness': 0.5, + 'erodeconvexdefects': 0, + 'minparticlearea': 20, + 'minparticledistance': 20, + 'measurefrac': 1, + 'compactness': 0.1, + 'seedRad': 3} + + self.ramanpoints = [] + self.particlecontours = [] + self.particlestats = [] + self.ramanscansortindex = None + self.ramanscandone = False + + self.results = {'polymers': None, + 'hqis': None, + 'additives': None, + 'additive_hqis': None} + + self.resultParams = {'minHQI': None, + 'compHQI': None} + self.spectraPath = None + self.particles2spectra = None #links idParticle to corresponding idSpectra (i.e., first measured particle (ID=0) is linked to spectra indices 0 and 1) + self.colorSeed = 'default' + self.resultsUploadedToSQL = [] + + self.readin = True # a value that is always set to True at loadData + # and mark that the coordinate system might be changed in the meantime + self.mode = "prepare" + if newProject: + self.fname = self.newProject(fname) + self.updatePath() + + def __eq__(self, other): + return recursiveDictCompare(self.__dict__, other.__dict__) + + def getPixelScale(self, mode=None): + if mode is None: + mode = self.imagescanMode + return (self.pixelscale_df if mode == 'df' else self.pixelscale_bf) + + def saveZvalImg(self): + if self.zvalimg is not None: + cv2imwrite_fix(self.getZvalImageName(), self.zvalimg) + + def loadZvalImg(self): + if os.path.exists(self.getZvalImageName()): + self.zvalimg = cv2imread_fix(self.getZvalImageName(), cv2.IMREAD_GRAYSCALE) + + def legacyConversion(self, recreatefullimage=False): + if self.version==0: + print("Converting legacy version 0 to 1") + print("This may take some time") + + # local imports as these functions are only needed for the rare occasion of legacy conversion + from opticalscan import loadAndPasteImage + + # try to load png and check for detection contours + recreatefullimage = recreatefullimage or not os.path.exists(self.getLegacyImageName()) + if not recreatefullimage: + img = cv2imread_fix(self.getLegacyImageName()) + Nc = len(self.particlecontours) + if Nc>0: + contour = self.particlecontours[Nc//2] + contpixels = img[contour[:,0,1],contour[:,0,0]] + if np.all(contpixels[:,1]==255) and np.all(contpixels[:,2]==0) \ + and np.all(contpixels[:,0]==0): + recreatefullimage = True + if not recreatefullimage: + cv2imwrite_fix(self.getImageName(), img) + del img + + if recreatefullimage: + print("recreating fullimage from grid data") + imgdata = None + zvalimg = None + Ngrid = len(self.grid) + + width, height, rotationvalue = self.imagedim_df + p0, p1 = self.maxdim[:2], self.maxdim[2:] + for i in range(Ngrid): + print(f"Processing image {i+1} of {Ngrid}") + names = [] + for k in range(len(self.zpositions)): + names.append(os.path.join(self.getScanPath(), f"image_{i}_{k}.bmp")) + p = self.grid[i] + imgdata, zvalimg = loadAndPasteImage(names, imgdata, zvalimg, width, + height, rotationvalue, p0, p1, p) + self.zvalimg = zvalimg + cv2imwrite_fix(self.getImageName(), cv2.cvtColor(imgdata, cv2.COLOR_RGB2BGR)) + del imgdata + self.saveZvalImg() + if "particleimgs" in self.__dict__: + del self.particleimgs + + self.version = 1 + + + if self.version == 1: + print("Converting legacy version 1 to 2") + if hasattr(self, 'pixelscale'): + print('pixelscale was', self.pixelscale) + self.pixelscale_bf = self.pixelscale + self.pixelscale_df = self.pixelscale + del self.pixelscale + + if hasattr(self, 'imagedim'): + self.imagedim_bf = self.imagedim + self.imagedim_df = self.imagedim + del self.imagedim + + self.version = 2 + + # add later conversion for higher version numbers here + + def getSubImage(self, img, index, draw=True): + contour = self.particlecontours[index] + x0, x1 = contour[:,0,0].min(), contour[:,0,0].max() + y0, y1 = contour[:,0,1].min(), contour[:,0,1].max() + subimg = img[y0:y1+1,x0:x1+1].copy() + if draw: + cv2.drawContours(subimg, [contour], -1, (0,255,0), 1) + return subimg + + def getZval(self, pixelpos): + assert self.zvalimg is not None + i, j = int(round(pixelpos[1])), int(round(pixelpos[0])) + if i>=self.zvalimg.shape[0]: + print('error in getZval:', self.zvalimg.shape, i, j) + i = self.zvalimg.shape[0]-1 + if j>=self.zvalimg.shape[1]: + print('error in getZval:', self.zvalimg.shape, i, j) + j = self.zvalimg.shape[1]-1 + zp = self.zvalimg[i,j] + z0, z1 = self.zpositions.min(), self.zpositions.max() + return zp/255.*(z1-z0) + z0 + + def mapHeight(self, x, y): + assert not self.readin + assert self.heightmap is not None + return self.heightmap[0]*x + self.heightmap[1]*y + self.heightmap[2] + + def mapToPixel(self, p, mode='df', force=False): + if not force: + assert not self.readin + p0 = copy(self.lastpos) + + if self.coordinatetransform is not None: + z = 0. if len(p)<3 else p[2] + T, pc = self.coordinatetransform + p = (np.dot(np.array([p[0], p[1], z])-pc, T.T)) + + if mode == 'df': + p0[0] -= self.signx*self.imagedim_df[0]/2 + p0[1] -= self.signy*self.imagedim_df[1]/2 + x, y = self.signx*(p[0] - p0[0])/self.pixelscale_df, self.signy*(p[1] - p0[1])/self.pixelscale_df + + elif mode == 'bf': + p0[0] -= self.signx*self.imagedim_bf[0]/2 + p0[1] += self.signy*self.imagedim_bf[1]/2 + x, y = self.signx*(p[0] - p0[0])/self.pixelscale_bf, self.signy*(p[1] - p0[1])/self.pixelscale_bf + + else: + raise ValueError(f'mapToPixel mode: {mode} not understood') + + return x, y + + def mapToLength(self, pixelpos, mode='df', force=False, returnz=False): + if not force: + assert not self.readin + p0 = copy(self.lastpos) + p0[0] += self.coordOffset[0] + p0[1] += self.coordOffset[1] + if mode == 'df': + p0[0] -= self.signx*self.imagedim_df[0]/2 + p0[1] -= self.signy*self.imagedim_df[1]/2 + x, y = (self.signx*pixelpos[0]*self.pixelscale_df + p0[0]), (p0[1] + self.signy*pixelpos[1]*self.pixelscale_df) + elif mode == 'bf': + p0[0] -= self.signx*self.imagedim_bf[0]/2 + p0[1] -= self.signy*self.imagedim_bf[1]/2 + x, y = (self.signx*pixelpos[0]*self.pixelscale_bf + p0[0]), (p0[1] + self.signy*pixelpos[1]*self.pixelscale_bf) + else: + raise ValueError(f'mapToLength mode: {mode} not understood') + + z = None + if (returnz and self.zvalimg is not None) or self.coordinatetransform is not None: + z = self.mapHeight(x, y) + z += self.getZval(pixelpos) + + if self.coordinatetransform is not None: + T, pc = self.coordinatetransform + x, y, z = (np.dot(np.array([x,y,z]), T) + pc) + + if returnz: + return x, y, z + return x, y + + def mapToLengthRaman(self, pixelpos, microscopeMode='df', noz=False): + p0x, p0y, z = self.mapToLength(pixelpos, mode=microscopeMode, returnz=True) + x, y = p0x + self.pshift[0], p0y + self.pshift[1] + return x, y, z + + def newProject(self, fname): + path = os.path.split(fname)[0] + name = os.path.splitext(os.path.basename(fname))[0] + newpath = os.path.join(path, name) + fname = os.path.join(newpath, name + ".pkl") + if not os.path.exists(newpath): + os.mkdir(newpath) # for new projects a directory will be created + elif os.path.exists(fname): # if this project is already there, load it instead + self.__dict__.update(loadData(fname).__dict__) + return fname + + def getScanPath(self): + scandir = os.path.join(self.path, "scanimages") + if not os.path.exists(scandir): + os.mkdir(scandir) + return scandir + + def updatePath(self): + self.path = os.path.split(self.fname)[0] + self.name = os.path.splitext(os.path.basename(self.fname))[0] + + def getImageName(self): + return os.path.join(self.path, 'fullimage.tif') + + def getZvalImageName(self): + return os.path.join(self.path, "zvalues.tif") + + def getLegacyImageName(self): + return os.path.join(self.path, "fullimage.png") + + def getLegacyDetectImageName(self): + return os.path.join(self.path, "detectimage.png") + + def getDetectImageName(self): + raise NotImplementedError("No longer implemented due to change in API") + + def getTmpImageName(self): + return os.path.join(self.path, "tmp.bmp") + + def saveParticleData(self): + print('Not saving ParticleData into text file...:\nThe current output format might be wrong, if multiple spectra per particle are present...') +# if len(self.ramanscansortindex)>0: +# data = [] +# pixelscale = (self.pixelscale_df if self.imagescanMode == 'df' else self.pixelscale_bf) +# for i in self.ramanscansortindex: +# data.append(list(self.ramanpoints[i])+list(self.particlestats[i])) +# data = np.array(data) +# data[:,0], data[:,1], z = self.mapToLengthRaman((data[:,0], data[:,1]), microscopeMode=self.imagescanMode, noz=True) +# data[:,2:7] *= pixelscale +# header = "x [µm], y [µm], length [µm], height [µm], length_ellipse [µm], height_ellipse [µm]" +# if data.shape[1]>6: +# header = header + ", area [µm^2]" +# data[:,6] *= pixelscale +# np.savetxt(os.path.join(self.path, "particledata.txt"), data, +# header=header) + + def save(self): + saveData(self, self.fname) + + def saveBackup(self): + inc = 0 + while True: + directory = os.path.dirname(self.fname) + filename = self.name + '_backup_' + str(inc) + '.pkl' + path = os.path.join(directory, filename) + if os.path.exists(path): + inc += 1 + else: + saveData(self, path) + return filename + diff --git a/helperfunctions.py b/helperfunctions.py index 38bffa4c1578e3893a903f1a9934d289356d10e0..3447a5a842b6277359b3db1c894dc674c7c16fb0 100644 --- a/helperfunctions.py +++ b/helperfunctions.py @@ -23,7 +23,16 @@ import numpy as np import cv2 import os +try: + from skimage.io import imread as skimread + from skimage.io import imsave as skimsave +except ImportError: + skimread = None + skimsave = None + def cv2imread_fix(fname, flags=cv2.IMREAD_COLOR): + if skimread is not None: + return skimread(fname, as_gray=(flags==cv2.IMREAD_GRAYSCALE)) with open(fname, "rb") as fp: cont = fp.read() img = cv2.imdecode(np.fromstring(cont, dtype=np.uint8), flags) @@ -31,6 +40,8 @@ def cv2imread_fix(fname, flags=cv2.IMREAD_COLOR): return None def cv2imwrite_fix(fname, img, params=None): + if skimsave is not None: + skimsave(fname, img) pathname, ext = os.path.splitext(fname) if params is None: ret, data = cv2.imencode(ext, img) diff --git a/zeissimporter.py b/zeissimporter.py index 8f2f54feb1e2df928615e14877b26f38f3580d96..c8c2c836f7b31f18db2b148b8526117cd9a7c787 100644 --- a/zeissimporter.py +++ b/zeissimporter.py @@ -1,270 +1,263 @@ -# -*- coding: utf-8 -*- -""" -GEPARD - Gepard-Enabled PARticle Detection -Copyright (C) 2018 Lars Bittrich and Josef Brandt, Leibniz-Institut für -Polymerforschung Dresden e. V. - -This program is free software: you can redistribute it and/or modify -it under the terms of the GNU General Public License as published by -the Free Software Foundation, either version 3 of the License, or -(at your option) any later version. - -This program is distributed in the hope that it will be useful, -but WITHOUT ANY WARRANTY; without even the implied warranty of -MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -GNU General Public License for more details. - -You should have received a copy of the GNU General Public License -along with this program, see COPYING. -If not, see . -""" -import os -from PyQt5 import QtCore, QtWidgets +# -*- coding: utf-8 -*- +""" +GEPARD - Gepard-Enabled PARticle Detection +Copyright (C) 2018 Lars Bittrich and Josef Brandt, Leibniz-Institut für +Polymerforschung Dresden e. V. + +This program is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 3 of the License, or +(at your option) any later version. + +This program is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with this program, see COPYING. +If not, see . +""" +import os +from PyQt5 import QtCore, QtWidgets from .zeissxml import ZeissHandler, make_parser from .opticalscan import PointCoordinates from .helperfunctions import cv2imread_fix, cv2imwrite_fix from .ramancom.ramancontrol import defaultPath from .dataset import DataSet -from scipy.optimize import least_squares -from itertools import permutations -import cv2 -import numpy as np - -class ZeissImporter(QtWidgets.QDialog): - def __init__(self, fname, ramanctrl, parent=None): - super().__init__(parent) - - if not ramanctrl.connect() or not self.readImportData(fname): - msg = QtWidgets.QMessageBox() - msg.setText('Connection failed! Please enable remote control.') - msg.exec() - self.validimport = False - return - else: - self.validimport = True - - self.ramanctrl = ramanctrl - - vbox = QtWidgets.QVBoxLayout() - pointgroup = QtWidgets.QGroupBox('Marker coordinates at Raman spot [µm]', self) - self.points = PointCoordinates(len(self.markers), self.ramanctrl, self, - names = [m.name for m in self.markers]) - self.points.pimageOnly.setVisible(False) - pointgroup.setLayout(self.points) - self.points.readPoint.connect(self.takePoint) - - self.pconvert = QtWidgets.QPushButton('Convert', self) - self.pexit = QtWidgets.QPushButton('Cancel', self) - self.pconvert.released.connect(self.convert) - self.pexit.released.connect(self.reject) - self.pconvert.setEnabled(False) - - btnLayout = QtWidgets.QHBoxLayout() - btnLayout.addStretch() - btnLayout.addWidget(self.pconvert) - btnLayout.addWidget(self.pexit) - - label = QtWidgets.QLabel("Z-Image blur radius", self) - self.blurspinbox = QtWidgets.QSpinBox(self) - self.blurspinbox.setMinimum(3) - self.blurspinbox.setMaximum(99) - self.blurspinbox.setSingleStep(2) - self.blurspinbox.setValue(5) - blurlayout = QtWidgets.QHBoxLayout() - blurlayout.addWidget(label) - blurlayout.addWidget(self.blurspinbox) - blurlayout.addStretch() - - vbox.addWidget(pointgroup) - vbox.addLayout(blurlayout) - vbox.addLayout(btnLayout) - self.setLayout(vbox) - - def readImportData(self, fname): - path = os.path.split(fname)[0] - self.zmapimgname = os.path.join(path, '3D.tif') - self.edfimgname = os.path.join(path, 'EDF.tif') - xmlname = os.path.join(path, '3D.tif_metadata.xml') - - errmsges = [] - if not os.path.exists(self.zmapimgname): - errmsges.append('Depth map image not found: 3D.tif') - if not os.path.exists(self.edfimgname): - errmsges.append('EDF image not found: EDF.tif') - if not os.path.exists(xmlname): - errmsges.append('XML metadata not found: 3D.tif_metadata.xml') - else: - parser = make_parser() - z = ZeissHandler() - parser.setContentHandler(z) - parser.parse(xmlname) - - if len(z.markers)<3: - errmsges.append('Fewer than 3 markers found to adjust coordinates!') - if None in [z.region.centerx, z.region.centery, - z.region.width, z.region.height]: - errmsges.append('Image dimensions incomplete or missing!') - if None in [z.zrange.z0, z.zrange.zn, z.zrange.dz]: - errmsges.append('ZStack information missing or incomplete!') - - if len(errmsges)>0: - QtWidgets.QMessageBox.error(self, 'Error!', - '\n'.join(errmsges), - QtWidgets.QMessageBox.Ok, - QtWidgets.QMessageBox.Ok) - return False - - self.region = z.region - self.zrange = z.zrange - self.markers = z.markers - return True - - @QtCore.pyqtSlot(float, float, float) - def takePoint(self, x, y, z): - points = self.points.getPoints() - if len(points)>=3: - self.pconvert.setEnabled(True) - - @QtCore.pyqtSlot() - def convert(self): - fname = QtWidgets.QFileDialog.getSaveFileName(self, - 'Create New GEPARD Project', defaultPath, '*.pkl')[0] - if fname=='': - return - dataset = DataSet(fname, newProject=True) - T, pc, zpc = self.getTransform() - imgshape, warp_mat = self.convertZimg(dataset, T, pc, zpc) - self.convertImage(dataset, warp_mat) - dataset.save() - self.gepardname = dataset.fname - self.accept() - - def convertImage(self, dataset, warp_mat): - img = cv2imread_fix(self.edfimgname) - img = cv2.warpAffine(img, warp_mat, img.shape[:2][::-1]) - cv2imwrite_fix(dataset.getImageName(), img) - - - def convertZimg(self, dataset, T, pc, zpc): - N = int(round((self.zrange.zn-self.zrange.z0)/self.zrange.dz)) - dataset.zpositions = np.linspace(self.zrange.z0, - self.zrange.zn, N)-zpc[2]+pc[2] - zimg = cv2imread_fix(self.zmapimgname, cv2.IMREAD_GRAYSCALE) - zmdist = zimg.mean() - zm = zmdist/255.*(self.zrange.zn-self.zrange.z0) + self.zrange.z0 - - radius = self.blurspinbox.value() - blur = cv2.GaussianBlur(zimg, (radius, radius), 0) - - pshift = self.ramanctrl.getRamanPositionShift() - dataset.pshift = pshift - pixelscale = self.region.width/zimg.shape[1] - - # use input image as single image aquired in one shot - dataset.imagedim_df = (self.region.width, self.region.height, 0.0) - dataset.pixelscale_df = pixelscale - - dataset.imagedim_bf = (self.region.width, self.region.height, 0.0) - dataset.pixelscale_bf = pixelscale - - # set image center as reference point in data set (transform from Zeiss) - p0 = np.dot((np.array([self.region.centerx, - self.region.centery,zm])-zpc),T)[:2] + pc[:2] - dataset.readin = False - dataset.lastpos = p0 - dataset.maxdim = p0 + p0 - - # pixel triangle for coordinate warping transformation - srcTri = np.array( [[0, 0], [zimg.shape[1] - 1, 0], - [0, zimg.shape[0] - 1]] ).astype(np.float32) - # upper left point (0,0) in Zeiss coordinates: - z0 = np.array([self.region.centerx - self.region.width/2, - self.region.centery + self.region.height/2]) - # transform pixel data to Zeiss coordinates - dstTri = np.array([[p[0]*pixelscale + z0[0], - z0[1] - p[1]*pixelscale, zm] for p in srcTri]).astype(np.double)-zpc - - # transform to Raman coordinates - dstTri = np.dot(dstTri,T) + pc[np.newaxis,:] - - # tilt blur image based on transformend z and adapt zpositions - x = np.linspace(0,1,blur.shape[1]) - y = np.linspace(0,1,blur.shape[0]) - x, y = np.meshgrid(x,y) - zmap = x*(dstTri[1,2]-dstTri[0,2]) + y*(dstTri[2,2]-dstTri[0,2]) + \ - (zimg * ((self.zrange.zn-self.zrange.z0)/255.) - \ - zmdist*((self.zrange.zn-self.zrange.z0)/255.)) - zmin, zmax = zmap.min(), zmap.max() - dataset.zpositions = np.array([zmap.min(), zmap.max()]) - blur = (zmap-zmin)*(255./(zmax-zmin)) - blur[blur>255.] = 255. - blur = np.uint8(blur) - # transform triangle back to pixel - dstTri = np.array([dataset.mapToPixel(p[:2]) for p in dstTri]).astype(np.float32) - - warp_mat = cv2.getAffineTransform(srcTri, dstTri) - - blur = cv2.warpAffine(blur, warp_mat, zimg.shape[::-1]) - zimgname = dataset.getZvalImageName() - cv2imwrite_fix(zimgname, blur) - return zimg.shape, warp_mat - - def getTransform(self): - points = self.points.getPoints() - pshift = self.ramanctrl.getRamanPositionShift() - points[:,0] -= pshift[0] - points[:,1] -= pshift[1] - zpoints = np.array([m.getPos() for m in self.markers], dtype=np.double) - pc = points.mean(axis=0) - zpc = zpoints.mean(axis=0) - - points -= pc[np.newaxis,:] - zpoints -= zpc[np.newaxis,:] - - def getRotMat(angles): - c1, s1 = np.cos(angles[0]), np.sin(angles[0]) - c2, s2 = np.cos(angles[1]), np.sin(angles[1]) - c3, s3 = np.cos(angles[2]), np.sin(angles[2]) - return np.mat([[c1*c3-s1*c2*s3, -c1*s3-s1*c2*c3, s1*s2], - [s1*c3+c1*c2*s3, -s1*s3+c1*c2*c3, -c1*s2], - [s1*s3, s2*c3, c2]]) - - # find the transformation matrix with best fit for small angles in - # [-45°,45°] for all permutation of markers - permbest = None - pointsbest = None - for perm in permutations(range(points.shape[0])): - ppoints = points[perm,:] - - def err(angles_shift): - T = getRotMat(angles_shift[:3]).T.A - return (np.dot(zpoints, T) - angles_shift[np.newaxis,3:] \ - - ppoints).ravel() - - angle = np.zeros(3) - opt = least_squares(err, np.concatenate((angle, np.zeros(3))), - bounds=(np.array([-np.pi/4]*3+[-np.inf]*3), - np.array([np.pi/4]*3+[np.inf]*3)), - method='dogbox') - - if permbest is None or \ - permbest.cost>opt.cost: - print("Current best permutation:", perm, flush=True) - permbest = opt - pointsbest = ppoints - - optangles = permbest.x[:3] - shift = permbest.x[3:] - T = getRotMat(optangles).T.A - e = (np.dot(zpoints, T)-shift[np.newaxis,:]-pointsbest) - print("Transformation angles:", optangles, flush=True) - print("Transformation shift:", shift, flush=True) - print("Transformation err:", e, flush=True) - d = np.linalg.norm(e, axis=1) - if np.any(d>1.): - QtWidgets.QMessageBox.warning(self, 'Warning!', - f'Transformation residuals are large:{d}', - QtWidgets.QMessageBox.Ok, - QtWidgets.QMessageBox.Ok) - return T, pc-shift, zpc - \ No newline at end of file +from scipy.optimize import least_squares +from itertools import permutations +import cv2 +import numpy as np + +class ZeissImporter(QtWidgets.QDialog): + def __init__(self, fname, ramanctrl, parent=None): + super().__init__(parent) + + if not ramanctrl.connect() or not self.readImportData(fname): + msg = QtWidgets.QMessageBox() + msg.setText('Connection failed! Please enable remote control.') + msg.exec() + self.validimport = False + return + else: + self.validimport = True + + self.ramanctrl = ramanctrl + + vbox = QtWidgets.QVBoxLayout() + pointgroup = QtWidgets.QGroupBox('Marker coordinates at Raman spot [µm]', self) + self.points = PointCoordinates(len(self.markers), self.ramanctrl, self, + names = [m.name for m in self.markers]) + self.points.pimageOnly.setVisible(False) + pointgroup.setLayout(self.points) + self.points.readPoint.connect(self.takePoint) + + self.pconvert = QtWidgets.QPushButton('Convert', self) + self.pexit = QtWidgets.QPushButton('Cancel', self) + self.pconvert.released.connect(self.convert) + self.pexit.released.connect(self.reject) + self.pconvert.setEnabled(False) + + self.xinvert = QtWidgets.QCheckBox('Invert x-axis') + self.yinvert = QtWidgets.QCheckBox('Invert y-axis') + self.zinvert = QtWidgets.QCheckBox('Invert z-axis') + + btnLayout = QtWidgets.QHBoxLayout() + btnLayout.addStretch() + btnLayout.addWidget(self.pconvert) + btnLayout.addWidget(self.pexit) + + label = QtWidgets.QLabel("Z-Image blur radius", self) + self.blurspinbox = QtWidgets.QSpinBox(self) + self.blurspinbox.setMinimum(3) + self.blurspinbox.setMaximum(99) + self.blurspinbox.setSingleStep(2) + self.blurspinbox.setValue(5) + blurlayout = QtWidgets.QHBoxLayout() + blurlayout.addWidget(label) + blurlayout.addWidget(self.blurspinbox) + blurlayout.addStretch() + + vbox.addWidget(pointgroup) + vbox.addLayout(blurlayout) + vbox.addWidget(self.xinvert) + vbox.addWidget(self.yinvert) + vbox.addWidget(self.zinvert) + vbox.addLayout(btnLayout) + self.setLayout(vbox) + + def readImportData(self, fname): + path = os.path.split(fname)[0] + self.edfimgname, self.zmapimgname, xmlname = '', '', '' + for name in os.listdir(path): + if name.lower().endswith('_meta.xml'): + xmlname = os.path.join(path, name) + elif name.lower().endswith('_c1.tif'): + self.edfimgname = os.path.join(path, name) + elif name.lower().endswith('_c2.tif'): + self.zmapimgname = os.path.join(path, name) + + errmsges = [] + if not os.path.exists(self.zmapimgname): + errmsges.append('Depth map image not found: NAME_c2.tif') + if not os.path.exists(self.edfimgname): + errmsges.append('EDF image not found: NAME_c1.tif') + if not os.path.exists(xmlname): + errmsges.append('XML metadata not found: NAME_meta.xml') + else: + parser = make_parser() + z = ZeissHandler() + parser.setContentHandler(z) + parser.parse(xmlname) + + if len(z.markers)<3: + errmsges.append('Fewer than 3 markers found to adjust coordinates!') + if None in [z.region.centerx, z.region.centery, + z.region.width, z.region.height]: + errmsges.append('Image dimensions incomplete or missing!') + if None in [z.zrange.z0, z.zrange.zn, z.zrange.dz]: + errmsges.append('ZStack information missing or incomplete!') + + if len(errmsges)>0: + QtWidgets.QMessageBox.critical(self, 'Error!', + '\n'.join(errmsges), + QtWidgets.QMessageBox.Ok, + QtWidgets.QMessageBox.Ok) + return False + + self.region = z.region + self.zrange = z.zrange + self.markers = z.markers + print(self.region) + print(self.zrange) + print(self.markers, flush=True) + return True + + @QtCore.pyqtSlot(float, float, float) + def takePoint(self, x, y, z): + points = self.points.getPoints() + if len(points)>=3: + self.pconvert.setEnabled(True) + + @QtCore.pyqtSlot() + def convert(self): + T, pc, zpc, accept = self.getTransform() + if accept: + fname = QtWidgets.QFileDialog.getSaveFileName(self, + 'Create New GEPARD Project', defaultPath, '*.pkl')[0] + if fname=='': + return + dataset = DataSet(fname, newProject=True) + self.convertZimg(dataset, T, pc, zpc) + self.convertImage(dataset) + dataset.save() + self.gepardname = dataset.fname + self.accept() + + def convertImage(self, dataset): + img = cv2imread_fix(self.edfimgname) + cv2imwrite_fix(dataset.getImageName(), img) + + def convertZimg(self, dataset, T, pc, zpc): + N = int(round(abs(self.zrange.zn-self.zrange.z0)/self.zrange.dz)) + z0, zn = self.zrange.z0, self.zrange.zn + if zn255.] = 255. + blur[np.isnan(blur)] = 0. + blur = np.uint8(blur) + + zimgname = dataset.getZvalImageName() + cv2imwrite_fix(zimgname, blur) + dataset.zvalimg = "saved" + + def getTransform(self): + points = self.points.getPoints() + pshift = self.ramanctrl.getRamanPositionShift() + points[:,0] -= pshift[0] + points[:,1] -= pshift[1] + Parity = np.mat(np.diag([-1. if self.xinvert.isChecked() else 1., + -1. if self.yinvert.isChecked() else 1., + -1. if self.zinvert.isChecked() else 1.])) + zpoints = np.array([m.getPos() for m in self.markers], dtype=np.double) + pc = points.mean(axis=0) + zpc = zpoints.mean(axis=0) + + points -= pc[np.newaxis,:] + zpoints -= zpc[np.newaxis,:] + + def getRotMat(angles): + c1, s1 = np.cos(angles[0]), np.sin(angles[0]) + c2, s2 = np.cos(angles[1]), np.sin(angles[1]) + c3, s3 = np.cos(angles[2]), np.sin(angles[2]) + return np.mat([[c1*c3-s1*c2*s3, -c1*s3-s1*c2*c3, s1*s2], + [s1*c3+c1*c2*s3, -s1*s3+c1*c2*c3, -c1*s2], + [s1*s3, s2*c3, c2]]) + + # find the transformation matrix with best fit for small angles in + # [-45°,45°] for all permutation of markers + permbest = None + pointsbest = None + ppoints = points[:,:].copy() + + def err(angles_shift): + T = (getRotMat(angles_shift[:3]).T*Parity).A + return (np.dot(zpoints, T) - angles_shift[np.newaxis,3:] \ + - ppoints).ravel() + + angle = np.zeros(3) + opt = least_squares(err, np.concatenate((angle, np.zeros(3))), + bounds=(np.array([-np.pi/4]*3+[-np.inf]*3), + np.array([np.pi/4]*3+[np.inf]*3)), + method='dogbox') + permbest = opt + pointsbest = ppoints + + optangles = permbest.x[:3] + shift = permbest.x[3:] + T = (getRotMat(optangles).T*Parity).A + e = (np.dot(zpoints, T)-shift[np.newaxis,:]-pointsbest) + print("Transformation angles:", optangles, flush=True) + print("Transformation shift:", shift, flush=True) + print("Transformation err:", e, flush=True) + d = np.linalg.norm(e, axis=1) + accept = True + if np.any(d>1.): + ret = QtWidgets.QMessageBox.warning(self, 'Warning!', + f'Transformation residuals are large:{d}', + QtWidgets.QMessageBox.Ok|QtWidgets.QMessageBox.Cancel, + QtWidgets.QMessageBox.Ok) + if ret==QtWidgets.QMessageBox.Cancel: + accept = False + return T, pc-shift, zpc, accept +