|
| 1 | +import configparser |
| 2 | +import typing |
| 3 | + |
| 4 | +from dataclasses import dataclass |
| 5 | +from numba import uint16, uint32, uint64, float32 |
| 6 | +from numba import njit, prange |
| 7 | +from numba.types import ListType, Array |
| 8 | + |
| 9 | +import numpy as np |
| 10 | + |
| 11 | +import slicer |
| 12 | + |
| 13 | + |
| 14 | +""" |
| 15 | +Normalization (just for reference) |
| 16 | +def norm(x, minimum, maximum): |
| 17 | + return minimum + (maximum - minimum) * x / 65535.0 |
| 18 | +
|
| 19 | +
|
| 20 | +def denorm(x, current_min, current_max, into_min=0, into_max=65535): |
| 21 | + return into_min + (into_max - into_min) * (x - current_min) / (current_max - current_min) |
| 22 | +""" |
| 23 | + |
| 24 | + |
| 25 | +@njit(parallel=True) |
| 26 | +def fastMapping( |
| 27 | + porosityImageArray: Array(float32, 3, "C"), |
| 28 | + labelsArray: Array(uint16, 3, "C"), |
| 29 | + maskArray: Array(uint16, 3, "C"), |
| 30 | + targets: ListType(uint16), |
| 31 | +): |
| 32 | + """ |
| 33 | + Same as: |
| 34 | + # porosityImageFinalArray = np.clip(porosityImageFloat, 0, 1, dtype=np.float32) |
| 35 | + # porosityImageFinalArray[self.__lastLabelsArray >= targets[1]] = np.float32(0) |
| 36 | + # porosityImageFinalArray[self.__lastLabelsArray == 0] = np.float32(0) |
| 37 | + # porosityImageFinalArray[self.__lastLabelsArray == targets[0]] = np.float32(1) |
| 38 | + # porosityImageFinalArray *= self.__lastMaskArray |
| 39 | + """ |
| 40 | + |
| 41 | + z, y, x = porosityImageArray.shape |
| 42 | + for i in prange(z): |
| 43 | + for j in range(y): |
| 44 | + for k in range(x): |
| 45 | + |
| 46 | + if maskArray[i, j, k] == 0: |
| 47 | + porosityImageArray[i, j, k] = np.float32(0) |
| 48 | + continue |
| 49 | + |
| 50 | + if labelsArray[i, j, k] == targets[0]: |
| 51 | + porosityImageArray[i, j, k] = np.float32(1) |
| 52 | + elif labelsArray[i, j, k] >= targets[1]: |
| 53 | + porosityImageArray[i, j, k] = np.float32(0) |
| 54 | + elif labelsArray[i, j, k] == 0: |
| 55 | + porosityImageArray[i, j, k] = np.float32(0) |
| 56 | + else: |
| 57 | + porosityImageArray[i, j, k] = max(0, min(1, porosityImageArray[i, j, k])) |
| 58 | + |
| 59 | + return porosityImageArray |
| 60 | + |
| 61 | + |
| 62 | +@njit |
| 63 | +def fast1DBinCountUInt64(array: Array(uint32, 1, "C")) -> Array(uint64, 1, "C"): |
| 64 | + max_ = np.max(array) |
| 65 | + out = np.zeros(max_ + 1, dtype=np.uint64) |
| 66 | + |
| 67 | + for i in range(len(array)): |
| 68 | + out[array[i]] += 1 |
| 69 | + |
| 70 | + return out |
| 71 | + |
| 72 | + |
| 73 | +def mapValue(values, x: typing.Union[float, typing.List[float]]): |
| 74 | + # Compute the attenuation factor for the modelled region |
| 75 | + accsum = np.cumsum(values) |
| 76 | + return np.interp(x, accsum / accsum[-1], np.arange(1, len(values) + 1)) |
| 77 | + |
| 78 | + |
| 79 | +@dataclass |
| 80 | +class DataVar: |
| 81 | + values: np.ndarray |
| 82 | + start: uint32 |
| 83 | + label: uint32 |
| 84 | + name: str |
| 85 | + threshold: np.ndarray = None |
| 86 | + color: tuple = (0, 0, 0) |
| 87 | + |
| 88 | + |
| 89 | +class SampleModel: |
| 90 | + def __init__(self, pcrRange=(0.0, 1.0)): |
| 91 | + self.variables: typing.Dict[str, DataVar] = {} |
| 92 | + self.pcrMin = np.float32(pcrRange[0]) |
| 93 | + self.pcrMax = np.float32(pcrRange[1]) |
| 94 | + self.image = None |
| 95 | + self.limits = None |
| 96 | + |
| 97 | + def addData(self, image: np.ndarray, mask: np.ndarray, label: int = 0, color: tuple = (0, 0, 0)): |
| 98 | + self.image = image |
| 99 | + self.addSubGroup("Data", mask, label, color, markers=[0.0001, 0.99]) |
| 100 | + self.limits = [int(np.round(v)) for v in self.variables["Data"].threshold] |
| 101 | + |
| 102 | + def addSubGroup(self, name: str, mask: np.ndarray, label: int, color: tuple, markers: typing.List[float] = None): |
| 103 | + if self.image is None: |
| 104 | + raise ValueError("Image is missing") |
| 105 | + |
| 106 | + if np.any(mask.shape != self.image.shape): |
| 107 | + raise ValueError("Mask and image must have the same shape") |
| 108 | + |
| 109 | + counts = fast1DBinCountUInt64(self.image[mask]) |
| 110 | + counts[0] = 0 # exclude background 0 |
| 111 | + threshold = mapValue(counts, x=markers) |
| 112 | + |
| 113 | + if self.limits is not None: |
| 114 | + min_, max_ = self.limits |
| 115 | + posmax = min(len(counts), max_) |
| 116 | + posmin = max(min_, 0) |
| 117 | + else: |
| 118 | + posmax = len(counts) |
| 119 | + posmin = 0 |
| 120 | + |
| 121 | + for i in range(posmin, posmax, 1): |
| 122 | + if counts[i] > 0: |
| 123 | + posmin = i |
| 124 | + break |
| 125 | + |
| 126 | + counts = counts if posmin == 0 else counts[posmin:] |
| 127 | + self.variables[name] = DataVar(counts, posmin, label, name, threshold, color) |
| 128 | + |
| 129 | + def threshold(self, key, normalized=False): |
| 130 | + value = self.variables[key].threshold |
| 131 | + if normalized: |
| 132 | + return self.pcrMin + (self.pcrMax - self.pcrMin) * value / 65535 |
| 133 | + return value |
| 134 | + |
| 135 | + def setThreshold(self, key, value: np.ndarray, normalized=False): |
| 136 | + if normalized: |
| 137 | + value = 65535 * (value - self.pcrMin) / (self.pcrMax - self.pcrMin) |
| 138 | + self.variables[key].threshold = value |
| 139 | + |
| 140 | + def getImage(self, normalized=False): |
| 141 | + if normalized: |
| 142 | + return self.pcrMin + (self.pcrMax - self.pcrMin) * self.image / 65535 |
| 143 | + return self.image |
0 commit comments