From 3f8b0d281b3058faf23dbd0363de440bd04c6574 Mon Sep 17 00:00:00 2001 From: Greg Burri Date: Tue, 26 Jan 2016 14:55:21 +0100 Subject: [PATCH] Split the module 'ImgTools' in many modules. --- Parasitemia/ParasitemiaCore/Classifier.fs | 3 +- Parasitemia/ParasitemiaCore/EEOver.fs | 7 +- Parasitemia/ParasitemiaCore/Ellipse.fs | 143 +------ Parasitemia/ParasitemiaCore/Granulometry.fs | 2 +- .../ParasitemiaCore/ImgTools/Drawing.fs | 79 ++++ Parasitemia/ParasitemiaCore/ImgTools/Edges.fs | 121 ++++++ .../ParasitemiaCore/ImgTools/Histogram.fs | 84 +++++ Parasitemia/ParasitemiaCore/ImgTools/IO.fs | 15 + .../ParasitemiaCore/ImgTools/ImgTools.fs | 49 +++ .../{ImgTools.fs => ImgTools/Morpho.fs} | 350 +----------------- Parasitemia/ParasitemiaCore/ImgTools/Otsu.fs | 47 +++ Parasitemia/ParasitemiaCore/MainAnalysis.fs | 47 +-- .../ParasitemiaCore/ParasitemiaCore.fsproj | 8 +- .../ParasitemiaCore/ParasitesMarker.fs | 3 + Parasitemia/ParasitemiaCore/UnitsOfMeasure.fs | 4 +- .../ParasitemiaUI/XAML/AboutWindow.xaml | 2 +- 16 files changed, 446 insertions(+), 518 deletions(-) create mode 100644 Parasitemia/ParasitemiaCore/ImgTools/Drawing.fs create mode 100644 Parasitemia/ParasitemiaCore/ImgTools/Edges.fs create mode 100644 Parasitemia/ParasitemiaCore/ImgTools/Histogram.fs create mode 100644 Parasitemia/ParasitemiaCore/ImgTools/IO.fs create mode 100644 Parasitemia/ParasitemiaCore/ImgTools/ImgTools.fs rename Parasitemia/ParasitemiaCore/{ImgTools.fs => ImgTools/Morpho.fs} (63%) create mode 100644 Parasitemia/ParasitemiaCore/ImgTools/Otsu.fs diff --git a/Parasitemia/ParasitemiaCore/Classifier.fs b/Parasitemia/ParasitemiaCore/Classifier.fs index 615840c..f4b0844 100644 --- a/Parasitemia/ParasitemiaCore/Classifier.fs +++ b/Parasitemia/ParasitemiaCore/Classifier.fs @@ -196,7 +196,6 @@ let findCells (ellipses: Ellipse list) (parasites: ParasitesMarker.Result) (img: then parasiteArea <- parasiteArea + 1 - let cellClass = if float darkStainPixels > config.Parameters.maxDarkStainRatio * (float nbElement) then @@ -204,7 +203,7 @@ let findCells (ellipses: Ellipse list) (parasites: ParasitesMarker.Result) (img: elif nucleusPixels.Count > 0 && parasiteArea >= minimumParasiteArea then - let infectionToRemove = ImgTools.connectedComponents parasites.parasite nucleusPixels + let infectionToRemove = Morpho.connectedComponents parasites.parasite nucleusPixels for p in infectionToRemove do nucleusData.[p.Y, p.X, 0] <- 0uy InfectedRBC diff --git a/Parasitemia/ParasitemiaCore/EEOver.fs b/Parasitemia/ParasitemiaCore/EEOver.fs index eb8ac0b..6862dc5 100644 --- a/Parasitemia/ParasitemiaCore/EEOver.fs +++ b/Parasitemia/ParasitemiaCore/EEOver.fs @@ -1,4 +1,5 @@ -module ParasitemiaCore.EEOver +// Translation from https://github.com/chraibi/EEOver. +module ParasitemiaCore.EEOver open System @@ -508,7 +509,9 @@ let private biquadroots (p: float[]) (r: float[,]) = quad () -// Return a tuple (area, x intersections, y intersections) +/// +/// Return a tuple (area, x intersections, y intersections). +/// let EEOverlapArea (e1: Types.Ellipse) (e2: Types.Ellipse) : (float32 * float32[] * float32[]) option = let h1, k1, a1, b1, phi_1 = float e1.Cx, float e1.Cy, float e1.A, float e1.B, float e1.Alpha let h2, k2, a2, b2, phi_2 = float e2.Cx, float e2.Cy, float e2.A, float e2.B, float e2.Alpha diff --git a/Parasitemia/ParasitemiaCore/Ellipse.fs b/Parasitemia/ParasitemiaCore/Ellipse.fs index 1059e3b..4e48230 100644 --- a/Parasitemia/ParasitemiaCore/Ellipse.fs +++ b/Parasitemia/ParasitemiaCore/Ellipse.fs @@ -14,142 +14,12 @@ open Config open MatchingEllipses open Const -type private SearchExtremum = Minimum | Maximum - -let private goldenSectionSearch (f: float -> float) (nbIter: int) (xmin: float) (xmax: float) (searchExtremum: SearchExtremum) : (float * float) = - let gr = 1. / 1.6180339887498948482 - let mutable a = xmin - let mutable b = xmax - let mutable c = b - gr * (b - a) - let mutable d = a + gr * (b - a) - - for i in 1 .. nbIter do - let mutable fc = f c - let mutable fd = f d - - if searchExtremum = Maximum - then - let tmp = fc - fc <- fd - fd <- tmp - - if fc < fd - then - b <- d - d <- c - c <- b - gr * (b - a) - else - a <- c - c <- d - d <- a + gr * (b - a) - - let x = (b + a) / 2. - x, f x - -// Ellipse.A is always equal or greater than Ellipse.B. -// Ellipse.Alpha is between 0 and Pi. +/// +/// Try to build an ellipse from three points and two tangents passing by the first and the second point. +/// 'Ellipse.A' is always equal or greater than Ellipse.B. +/// 'Ellipse.Alpha' is between 0 and Pi. +/// let ellipse (p1x: float) (p1y: float) (m1: float) (p2x: float) (p2y: float) (m2: float) (p3x: float) (p3y: float) : Types.Ellipse option = - let accuracy_extremum_search_1 = 10 // 3 - let accuracy_extremum_search_2 = 10 // 4 - - // p3 as the referencial. - let p1x = p1x - p3x - let p1y = p1y - p3y - - let p2x = p2x - p3x - let p2y = p2y - p3y - - // Convert to polar coordinates. - let alpha1 = atan m1 - let alpha2 = atan m2 - - let r1 = sqrt (p1x ** 2. + p1y ** 2.) - let theta1 = atan2 p1y p1x - - let r2 = sqrt (p2x ** 2. + p2y ** 2.) - let theta2 = atan2 p2y p2x - - let valid = - 4. * sin (alpha1 - theta1) * (-r1 * sin (alpha1 - theta1) + r2 * sin (alpha1 - theta2)) * - sin (alpha2 - theta2) * (-r1 * sin (alpha2 - theta1) + r2 * sin (alpha2 - theta2)) + - r1 * r2 * sin (alpha1 - alpha2) ** 2. * sin (theta1 - theta2) ** 2. < 0. - - if valid - then - let r theta = - (r1 * r2 * (r1 * (cos (alpha2 + theta - theta1 - theta2) - cos (alpha2 - theta) * cos (theta1 - theta2)) * sin (alpha1 - theta1) + r2 * (-cos (alpha1 + theta - theta1 - theta2) + cos (alpha1 - theta) * cos (theta1 - theta2)) * sin (alpha2 - theta2)) * sin (theta1 - theta2)) / - (sin (alpha1 - theta1) * sin (alpha2 - theta2) * (r1 * sin (theta - theta1) - r2 * sin (theta - theta2)) ** 2. - r1 * r2 * sin (alpha1 - theta) * sin (alpha2 - theta) * sin (theta1 - theta2) ** 2.) - - let rabs = r >> abs - - // We search for an interval [theta_a, theta_b] and assume the function is unimodal in this interval. - let thetaTan, _ = goldenSectionSearch rabs accuracy_extremum_search_1 0. Math.PI Maximum - let rTan = r thetaTan - - let PTanx = rTan * cos thetaTan - let PTany = rTan * sin thetaTan - - let d1a = tan alpha1 - let d1b = -d1a * p1x + p1y - - let d2a = tan alpha2 - let d2b = -d2a * p2x + p2y - - let d3a = -1. / tan thetaTan - let d3b = -d3a * PTanx + PTany - - let Ux = -(d1b - d2b) / (d1a - d2a) - let Uy = -(d2a * d1b - d1a * d2b) / (d1a - d2a) - - let Vx = -(d1b - d3b) / (d1a - d3a) - let Vy = -(d3a * d1b - d1a * d3b) / (d1a - d3a) - - let Wx = p1x + (p2x - p1x) / 2. - let Wy = p1y + (p2y - p1y) / 2. - - let Zx = p1x + (PTanx - p1x) / 2. - let Zy = p1y + (PTany - p1y) / 2. - - let va = -(-Vy + Zy) / (Vx - Zx) - let vb = -(Zx * Vy - Vx * Zy) / (Vx - Zx) - - let ua = -(-Uy + Wy) / (Ux - Wx) - let ub = -(Wx * Uy - Ux * Wy) / (Ux - Wx) - - let cx = -(vb - ub) / (va - ua) - let cy = -(ua * vb - va * ub) / (va - ua) - - let rc = sqrt (cx ** 2. + cy ** 2.) - let psi = atan2 cy cx - - let rellipse theta = - sqrt ( - rc ** 2. + (r1 ** 2. * r2 ** 2. * (r1 * (cos (alpha2 + theta - theta1 - theta2) - cos (alpha2 - theta) * cos (theta1 - theta2)) * sin (alpha1 - theta1) + r2 * (-cos (alpha1 + theta - theta1 - theta2) + cos (alpha1 - theta) * cos (theta1 - theta2)) * sin (alpha2 - theta2)) ** 2. * sin (theta1 - theta2) ** 2.) / - (sin (alpha1 - theta1) * sin (alpha2 - theta2) * (r1 * sin (theta - theta1) - r2 * sin (theta - theta2)) ** 2. - r1 * r2 * sin (alpha1 - theta) * sin (alpha2 - theta) * sin (theta1 - theta2) ** 2.) ** 2. - - (2. * r1 * r2 * rc * cos (theta - psi) * (r1 * (cos (alpha2 + theta - theta1 - theta2) - cos (alpha2 - theta) * cos (theta1 - theta2)) * sin (alpha1 - theta1) + r2 * (-cos (alpha1 + theta - theta1 - theta2) + cos (alpha1 - theta) * cos (theta1 - theta2)) * sin (alpha2 - theta2)) * sin (theta1 - theta2)) / - (sin (alpha1 - theta1) * sin (alpha2 - theta2) * (r1 * sin (theta - theta1) - r2 * sin (theta - theta2)) ** 2. - r1 * r2 * sin (alpha1 - theta) * sin (alpha2 - theta) * sin (theta1 - theta2) ** 2.)) - - // We search for an interval [theta_a, theta_b] and assume the function is unimodal in this interval. - let r1eTheta, r1e = goldenSectionSearch rellipse accuracy_extremum_search_2 0. (Math.PI / 2.) Maximum // Pi/2 and not pi because the period is Pi. - let r2eTheta, r2e = goldenSectionSearch rellipse accuracy_extremum_search_2 0. (Math.PI / 2.) Minimum - - let rr1e = r r1eTheta - let r1ex = rr1e * cos r1eTheta - let r1ey = rr1e * sin r1eTheta - let mutable alpha = atan ((r1ey - cy) / (r1ex - cx)) - if alpha < 0. - then - alpha <- alpha + Math.PI - - // Ride off the p3 referential. - let cx = cx + p3x - let cy = cy + p3y - - Some (Types.Ellipse(float32 cx, float32 cy, float32 r1e, float32 r2e, float32 alpha)) - else - None - -let ellipse2 (p1x: float) (p1y: float) (m1: float) (p2x: float) (p2y: float) (m2: float) (p3x: float) (p3y: float) : Types.Ellipse option = let p0 = pointFromTwoLines (Types.Line(float32 m1, float32 (p1y - m1 * p1x))) (Types.Line(float32 m2, float32(p2y - m2 * p2x))) let p0x, p0y = float p0.X, float p0.Y @@ -259,7 +129,6 @@ let private areVectorsValid (p1x: float32) (p1y: float32) (p2x: float32) (p2y: f else Some (m1, m2) - let find (edges: Matrix) (xGradient: Matrix) (yGradient: Matrix) @@ -336,7 +205,7 @@ let find (edges: Matrix) then match areVectorsValid (float32 p1xf) (float32 p1yf) (float32 p2xf) (float32 p2yf) -xDirData.[p1.Y, p1.X] -yDirData.[p1.Y, p1.X] -xDirData.[p2.Y, p2.X] -yDirData.[p2.Y, p2.X] with | Some (m1, m2) -> - match ellipse2 p1xf p1yf (float m1) p2xf p2yf (float m2) p3xf p3yf with + match ellipse p1xf p1yf (float m1) p2xf p2yf (float m2) p3xf p3yf with | Some e when e.Cx > 0.f && e.Cx < w_f - 1.f && e.Cy > 0.f && e.Cy < h_f - 1.f && e.A >= r1 - radiusTolerance && e.A <= r2 + radiusTolerance && e.B >= r1 - radiusTolerance && e.B <= r2 + radiusTolerance -> ellipses.Add e diff --git a/Parasitemia/ParasitemiaCore/Granulometry.fs b/Parasitemia/ParasitemiaCore/Granulometry.fs index 0afe5cc..b325534 100644 --- a/Parasitemia/ParasitemiaCore/Granulometry.fs +++ b/Parasitemia/ParasitemiaCore/Granulometry.fs @@ -71,7 +71,7 @@ let findRadiusByAreaClosing (img: Image) (radiusRange: int * int) let mutable maxDiff = 0.f let mutable max_r = r1 - ImgTools.areaCloseFWithFun imgCopy [ for r in r1 .. r2 -> Math.PI * float r ** 2. |> roundInt, r ] (fun r diff -> + Morpho.areaCloseFWithFun imgCopy [ for r in r1 .. r2 -> Math.PI * float r ** 2. |> roundInt, r ] (fun r diff -> if r <> r1 && diff > maxDiff then maxDiff <- diff diff --git a/Parasitemia/ParasitemiaCore/ImgTools/Drawing.fs b/Parasitemia/ParasitemiaCore/ImgTools/Drawing.fs new file mode 100644 index 0000000..363acda --- /dev/null +++ b/Parasitemia/ParasitemiaCore/ImgTools/Drawing.fs @@ -0,0 +1,79 @@ +module ParasitemiaCore.Drawing + +open System +open System.Drawing + +open Emgu.CV +open Emgu.CV.Structure + +open Const +open Types + +let drawPoints (img: Image) (points: Points) (intensity: 'TDepth) = + for p in points do + img.Data.[p.Y, p.X, 0] <- intensity + +let drawLine (img: Image<'TColor, 'TDepth>) (color: 'TColor) (x0: int) (y0: int) (x1: int) (y1: int) (thickness: int) = + img.Draw(LineSegment2D(Point(x0, y0), Point(x1, y1)), color, thickness); + +let drawLineF (img: Image<'TColor, 'TDepth>) (color: 'TColor) (x0: float) (y0: float) (x1: float) (y1: float) (thickness: int) = + img.Draw(LineSegment2DF(PointF(float32 x0, float32 y0), PointF(float32 x1, float32 y1)), color, thickness, CvEnum.LineType.AntiAlias); + +let drawEllipse (img: Image<'TColor, 'TDepth>) (e: Ellipse) (color: 'TColor) (alpha: float) = + if alpha >= 1.0 + then + img.Draw(Emgu.CV.Structure.Ellipse(PointF(e.Cx, e.Cy), SizeF(2.f * e.B, 2.f * e.A), e.Alpha / PI * 180.f), color, 1, CvEnum.LineType.AntiAlias) + else + let windowPosX = e.Cx - e.A - 5.f + let gapX = windowPosX - (float32 (int windowPosX)) + + let windowPosY = e.Cy - e.A - 5.f + let gapY = windowPosY - (float32 (int windowPosY)) + + let roi = Rectangle(int windowPosX, int windowPosY, 2.f * (e.A + 5.f) |> int, 2.f * (e.A + 5.f) |> int) + + img.ROI <- roi + if roi = img.ROI // We do not display ellipses touching the edges (FIXME) + then + use i = new Image<'TColor, 'TDepth>(img.ROI.Size) + i.Draw(Emgu.CV.Structure.Ellipse(PointF(e.A + 5.f + gapX, e.A + 5.f + gapY), SizeF(2.f * e.B, 2.f * e.A), e.Alpha / PI * 180.f), color, 1, CvEnum.LineType.AntiAlias) + CvInvoke.AddWeighted(img, 1.0, i, alpha, 0.0, img) + img.ROI <- Rectangle.Empty + +let drawEllipses (img: Image<'TColor, 'TDepth>) (ellipses: Ellipse list) (color: 'TColor) (alpha: float) = + List.iter (fun e -> drawEllipse img e color alpha) ellipses + +let rngCell = System.Random() +let drawCell (img: Image) (drawCellContent: bool) (c: Cell) = + if drawCellContent + then + let colorB = rngCell.Next(20, 70) + let colorG = rngCell.Next(20, 70) + let colorR = rngCell.Next(20, 70) + + for y in 0 .. c.elements.Height - 1 do + for x in 0 .. c.elements.Width - 1 do + if c.elements.[y, x] > 0uy + then + let dx, dy = c.center.X - c.elements.Width / 2, c.center.Y - c.elements.Height / 2 + let b = img.Data.[y + dy, x + dx, 0] |> int + let g = img.Data.[y + dy, x + dx, 1] |> int + let r = img.Data.[y + dy, x + dx, 2] |> int + img.Data.[y + dy, x + dx, 0] <- if b + colorB > 255 then 255uy else byte (b + colorB) + img.Data.[y + dy, x + dx, 1] <- if g + colorG > 255 then 255uy else byte (g + colorG) + img.Data.[y + dy, x + dx, 2] <- if r + colorR > 255 then 255uy else byte (r + colorR) + + let crossColor, crossColor2 = + match c.cellClass with + | HealthyRBC -> Bgr(255., 0., 0.), Bgr(255., 255., 255.) + | InfectedRBC -> Bgr(0., 0., 255.), Bgr(120., 120., 255.) + | Peculiar -> Bgr(0., 0., 0.), Bgr(80., 80., 80.) + + drawLine img crossColor2 (c.center.X - 3) c.center.Y (c.center.X + 3) c.center.Y 2 + drawLine img crossColor2 c.center.X (c.center.Y - 3) c.center.X (c.center.Y + 3) 2 + + drawLine img crossColor (c.center.X - 3) c.center.Y (c.center.X + 3) c.center.Y 1 + drawLine img crossColor c.center.X (c.center.Y - 3) c.center.X (c.center.Y + 3) 1 + +let drawCells (img: Image) (drawCellContent: bool) (cells: Cell list) = + List.iter (fun c -> drawCell img drawCellContent c) cells \ No newline at end of file diff --git a/Parasitemia/ParasitemiaCore/ImgTools/Edges.fs b/Parasitemia/ParasitemiaCore/ImgTools/Edges.fs new file mode 100644 index 0000000..b174ee9 --- /dev/null +++ b/Parasitemia/ParasitemiaCore/ImgTools/Edges.fs @@ -0,0 +1,121 @@ +module ParasitemiaCore.Edges + +open System +open System.Drawing +open System.Collections.Generic + +open Emgu.CV +open Emgu.CV.Structure + +open Const +open Histogram +open Otsu + +/// +/// Find edges of an image by using the Canny approach. +/// The thresholds are automatically defined with otsu on gradient magnitudes. +/// +/// +let find (img: Image) : Matrix * Matrix * Matrix = + let w = img.Width + let h = img.Height + + use sobelKernel = + new Matrix(array2D [[ 1.0f; 0.0f; -1.0f ] + [ 2.0f; 0.0f; -2.0f ] + [ 1.0f; 0.0f; -1.0f ]]) + + let xGradient = new Matrix(img.Size) + let yGradient = new Matrix(img.Size) + CvInvoke.Filter2D(img, xGradient, sobelKernel, Point(1, 1)) + CvInvoke.Filter2D(img, yGradient, sobelKernel.Transpose(), Point(1, 1)) + + use magnitudes = new Matrix(xGradient.Size) + use angles = new Matrix(xGradient.Size) + CvInvoke.CartToPolar(xGradient, yGradient, magnitudes, angles) // Compute the magnitudes and angles. + + let thresholdHigh, thresholdLow = + let sensibilityHigh = 0.1f + let sensibilityLow = 0.0f + let threshold, _, _ = otsu (histogramMat magnitudes 300) + threshold + (sensibilityHigh * threshold), threshold - (sensibilityLow * threshold) + + // Non-maximum suppression. + use nms = new Matrix(xGradient.Size) + + let nmsData = nms.Data + let anglesData = angles.Data + let magnitudesData = magnitudes.Data + let xGradientData = xGradient.Data + let yGradientData = yGradient.Data + + for i in 0 .. h - 1 do + nmsData.[i, 0] <- 0uy + nmsData.[i, w - 1] <- 0uy + + for j in 0 .. w - 1 do + nmsData.[0, j] <- 0uy + nmsData.[h - 1, j] <- 0uy + + for i in 1 .. h - 2 do + for j in 1 .. w - 2 do + let vx = xGradientData.[i, j] + let vy = yGradientData.[i, j] + if vx <> 0.f || vy <> 0.f + then + let angle = anglesData.[i, j] + + let vx', vy' = abs vx, abs vy + let ratio2 = if vx' > vy' then vy' / vx' else vx' / vy' + let ratio1 = 1.f - ratio2 + + let mNeigbors (sign: int) : float32 = + if angle < PI / 4.f + then ratio1 * magnitudesData.[i, j + sign] + ratio2 * magnitudesData.[i + sign, j + sign] + elif angle < PI / 2.f + then ratio2 * magnitudesData.[i + sign, j + sign] + ratio1 * magnitudesData.[i + sign, j] + elif angle < 3.f * PI / 4.f + then ratio1 * magnitudesData.[i + sign, j] + ratio2 * magnitudesData.[i + sign, j - sign] + elif angle < PI + then ratio2 * magnitudesData.[i + sign, j - sign] + ratio1 * magnitudesData.[i, j - sign] + elif angle < 5.f * PI / 4.f + then ratio1 * magnitudesData.[i, j - sign] + ratio2 * magnitudesData.[i - sign, j - sign] + elif angle < 3.f * PI / 2.f + then ratio2 * magnitudesData.[i - sign, j - sign] + ratio1 * magnitudesData.[i - sign, j] + elif angle < 7.f * PI / 4.f + then ratio1 * magnitudesData.[i - sign, j] + ratio2 * magnitudesData.[i - sign, j + sign] + else ratio2 * magnitudesData.[i - sign, j + sign] + ratio1 * magnitudesData.[i, j + sign] + + let m = magnitudesData.[i, j] + if m >= thresholdLow && m > mNeigbors 1 && m > mNeigbors -1 + then + nmsData.[i, j] <- 1uy + + // suppressMConnections nms // It's not helpful for the rest of the process (ellipse detection). + + let edges = new Matrix(xGradient.Size) + let edgesData = edges.Data + + // Hysteresis thresholding. + let toVisit = Stack() + for i in 0 .. h - 1 do + for j in 0 .. w - 1 do + if nmsData.[i, j] = 1uy && magnitudesData.[i, j] >= thresholdHigh + then + nmsData.[i, j] <- 0uy + toVisit.Push(Point(j, i)) + while toVisit.Count > 0 do + let p = toVisit.Pop() + edgesData.[p.Y, p.X] <- 1uy + for i' in -1 .. 1 do + for j' in -1 .. 1 do + if i' <> 0 || j' <> 0 + then + let ni = p.Y + i' + let nj = p.X + j' + if ni >= 0 && ni < h && nj >= 0 && nj < w && nmsData.[ni, nj] = 1uy + then + nmsData.[ni, nj] <- 0uy + toVisit.Push(Point(nj, ni)) + + edges, xGradient, yGradient \ No newline at end of file diff --git a/Parasitemia/ParasitemiaCore/ImgTools/Histogram.fs b/Parasitemia/ParasitemiaCore/ImgTools/Histogram.fs new file mode 100644 index 0000000..92b73f6 --- /dev/null +++ b/Parasitemia/ParasitemiaCore/ImgTools/Histogram.fs @@ -0,0 +1,84 @@ +module ParasitemiaCore.Histogram + +open System +open System.Drawing + +open Emgu.CV +open Emgu.CV.Structure + +type Histogram = { + data: int[] + total: int // Number of elements. + sum: int // Sum of all intensity. + min: float32 + max: float32 } + +let histogramImg (img: Image) (nbSamples: int) : Histogram = + let imgData = img.Data + + let min, max = + let min = ref [| 0.0 |] + let minLocation = ref <| [| Point() |] + let max = ref [| 0.0 |] + let maxLocation = ref <| [| Point() |] + img.MinMax(min, max, minLocation, maxLocation) + float32 (!min).[0], float32 (!max).[0] + + let inline bin (x: float32) : int = + let p = int ((x - min) / (max - min) * float32 nbSamples) + if p >= nbSamples then nbSamples - 1 else p + + let data = Array.zeroCreate nbSamples + + for i in 0 .. img.Height - 1 do + for j in 0 .. img.Width - 1 do + let p = bin imgData.[i, j, 0] + data.[p] <- data.[p] + 1 + + { data = data; total = img.Height * img.Width; sum = Array.sum data; min = min; max = max } + +let histogramMat (mat: Matrix) (nbSamples: int) : Histogram = + let matData = mat.Data + + let min, max = + let min = ref 0.0 + let minLocation = ref <| Point() + let max = ref 0.0 + let maxLocation = ref <| Point() + mat.MinMax(min, max, minLocation, maxLocation) + float32 !min, float32 !max + + let inline bin (x: float32) : int = + let p = int ((x - min) / (max - min) * float32 nbSamples) + if p >= nbSamples then nbSamples - 1 else p + + let data = Array.zeroCreate nbSamples + + for i in 0 .. mat.Height - 1 do + for j in 0 .. mat.Width - 1 do + let p = bin matData.[i, j] + data.[p] <- data.[p] + 1 + + { data = data; total = mat.Height * mat.Width; sum = Array.sum data; min = min; max = max } + +let histogram (values: float32 seq) (nbSamples: int) : Histogram = + let mutable min = Single.MaxValue + let mutable max = Single.MinValue + let mutable n = 0 + + for v in values do + n <- n + 1 + if v < min then min <- v + if v > max then max <- v + + let inline bin (x: float32) : int = + let p = int ((x - min) / (max - min) * float32 nbSamples) + if p >= nbSamples then nbSamples - 1 else p + + let data = Array.zeroCreate nbSamples + + for v in values do + let p = bin v + data.[p] <- data.[p] + 1 + + { data = data; total = n; sum = Array.sum data; min = min; max = max } \ No newline at end of file diff --git a/Parasitemia/ParasitemiaCore/ImgTools/IO.fs b/Parasitemia/ParasitemiaCore/ImgTools/IO.fs new file mode 100644 index 0000000..219eab6 --- /dev/null +++ b/Parasitemia/ParasitemiaCore/ImgTools/IO.fs @@ -0,0 +1,15 @@ +module ParasitemiaCore.IO + +open System +open System.Drawing + +open Emgu.CV +open Emgu.CV.Structure + +let saveImg (img: Image<'TColor, 'TDepth>) (filepath: string) = + img.Save(filepath) + +let saveMat (mat: Matrix<'TDepth>) (filepath: string) = + use img = new Image(mat.Size) + mat.CopyTo(img) + saveImg img filepath diff --git a/Parasitemia/ParasitemiaCore/ImgTools/ImgTools.fs b/Parasitemia/ParasitemiaCore/ImgTools/ImgTools.fs new file mode 100644 index 0000000..0510bc8 --- /dev/null +++ b/Parasitemia/ParasitemiaCore/ImgTools/ImgTools.fs @@ -0,0 +1,49 @@ +module ParasitemiaCore.ImgTools + +open System +open System.Drawing + +open Emgu.CV +open Emgu.CV.Structure + +let normalize (img: Image) (upperLimit: float) : Image = + let min = ref [| 0.0 |] + let minLocation = ref <| [| Point() |] + let max = ref [| 0.0 |] + let maxLocation = ref <| [| Point() |] + img.MinMax(min, max, minLocation, maxLocation) + let normalized = (img - (!min).[0]) / ((!max).[0] - (!min).[0]) + if upperLimit = 1.0 + then normalized + else upperLimit * normalized + +let mergeChannels (img: Image) (rgbWeights: float * float * float) : Image = + match rgbWeights with + | 1., 0., 0. -> img.[2] + | 0., 1., 0. -> img.[1] + | 0., 0., 1. -> img.[0] + | redFactor, greenFactor, blueFactor -> + let result = new Image(img.Size) + CvInvoke.AddWeighted(result, 1., img.[2], redFactor, 0., result) + CvInvoke.AddWeighted(result, 1., img.[1], greenFactor, 0., result) + CvInvoke.AddWeighted(result, 1., img.[0], blueFactor, 0., result) + result + +let mergeChannelsWithProjection (img: Image) (v1r: float32, v1g: float32, v1b: float32) (v2r: float32, v2g: float32, v2b: float32) (upperLimit: float) : Image = + let vr, vg, vb = v2r - v1r, v2g - v1g, v2b - v1b + let vMagnitude = sqrt (vr ** 2.f + vg ** 2.f + vb ** 2.f) + let project (r: float32) (g: float32) (b: float32) = ((r - v1r) * vr + (g - v1g) * vg + (b - v1b) * vb) / vMagnitude + let result = new Image(img.Size) + // TODO: Essayer en bindant Data pour gagner du temps + for i in 0 .. img.Height - 1 do + for j in 0 .. img.Width - 1 do + result.Data.[i, j, 0] <- project img.Data.[i, j, 2] img.Data.[i, j, 1] img.Data.[i, j, 0] + normalize result upperLimit + +// Normalize image values between 0uy and 255uy. +let normalizeAndConvert (img: Image) : Image = + (normalize (img.Convert()) 255.).Convert() + +let gaussianFilter (img : Image<'TColor, 'TDepth>) (standardDeviation : float) : Image<'TColor, 'TDepth> = + let size = 2 * int (ceil (4.0 * standardDeviation)) + 1 + img.SmoothGaussian(size, size, standardDeviation, standardDeviation) diff --git a/Parasitemia/ParasitemiaCore/ImgTools.fs b/Parasitemia/ParasitemiaCore/ImgTools/Morpho.fs similarity index 63% rename from Parasitemia/ParasitemiaCore/ImgTools.fs rename to Parasitemia/ParasitemiaCore/ImgTools/Morpho.fs index cbb65b4..b3a2e75 100644 --- a/Parasitemia/ParasitemiaCore/ImgTools.fs +++ b/Parasitemia/ParasitemiaCore/ImgTools/Morpho.fs @@ -1,4 +1,4 @@ -module ParasitemiaCore.ImgTools +module ParasitemiaCore.Morpho open System open System.Drawing @@ -8,172 +8,7 @@ open System.Linq open Emgu.CV open Emgu.CV.Structure -open Heap -open Const open Types -open Utils - -let normalize (img: Image) (upperLimit: float) : Image = - let min = ref [| 0.0 |] - let minLocation = ref <| [| Point() |] - let max = ref [| 0.0 |] - let maxLocation = ref <| [| Point() |] - img.MinMax(min, max, minLocation, maxLocation) - let normalized = (img - (!min).[0]) / ((!max).[0] - (!min).[0]) - if upperLimit = 1.0 - then normalized - else upperLimit * normalized - -let mergeChannels (img: Image) (rgbWeights: float * float * float) : Image = - match rgbWeights with - | 1., 0., 0. -> img.[2] - | 0., 1., 0. -> img.[1] - | 0., 0., 1. -> img.[0] - | redFactor, greenFactor, blueFactor -> - let result = new Image(img.Size) - CvInvoke.AddWeighted(result, 1., img.[2], redFactor, 0., result) - CvInvoke.AddWeighted(result, 1., img.[1], greenFactor, 0., result) - CvInvoke.AddWeighted(result, 1., img.[0], blueFactor, 0., result) - result - -let mergeChannelsWithProjection (img: Image) (v1r: float32, v1g: float32, v1b: float32) (v2r: float32, v2g: float32, v2b: float32) (upperLimit: float) : Image = - let vr, vg, vb = v2r - v1r, v2g - v1g, v2b - v1b - let vMagnitude = sqrt (vr ** 2.f + vg ** 2.f + vb ** 2.f) - let project (r: float32) (g: float32) (b: float32) = ((r - v1r) * vr + (g - v1g) * vg + (b - v1b) * vb) / vMagnitude - let result = new Image(img.Size) - // TODO: Essayer en bindant Data pour gagner du temps - for i in 0 .. img.Height - 1 do - for j in 0 .. img.Width - 1 do - result.Data.[i, j, 0] <- project img.Data.[i, j, 2] img.Data.[i, j, 1] img.Data.[i, j, 0] - normalize result upperLimit - -// Normalize image values between 0uy and 255uy. -let normalizeAndConvert (img: Image) : Image = - (normalize (img.Convert()) 255.).Convert() - -let saveImg (img: Image<'TColor, 'TDepth>) (filepath: string) = - img.Save(filepath) - -let saveMat (mat: Matrix<'TDepth>) (filepath: string) = - use img = new Image(mat.Size) - mat.CopyTo(img) - saveImg img filepath - -type Histogram = { data: int[]; total: int; sum: int; min: float32; max: float32 } - -let histogramImg (img: Image) (nbSamples: int) : Histogram = - let imgData = img.Data - - let min, max = - let min = ref [| 0.0 |] - let minLocation = ref <| [| Point() |] - let max = ref [| 0.0 |] - let maxLocation = ref <| [| Point() |] - img.MinMax(min, max, minLocation, maxLocation) - float32 (!min).[0], float32 (!max).[0] - - let inline bin (x: float32) : int = - let p = int ((x - min) / (max - min) * float32 nbSamples) - if p >= nbSamples then nbSamples - 1 else p - - let data = Array.zeroCreate nbSamples - - for i in 0 .. img.Height - 1 do - for j in 0 .. img.Width - 1 do - let p = bin imgData.[i, j, 0] - data.[p] <- data.[p] + 1 - - { data = data; total = img.Height * img.Width; sum = Array.sum data; min = min; max = max } - -let histogramMat (mat: Matrix) (nbSamples: int) : Histogram = - let matData = mat.Data - - let min, max = - let min = ref 0.0 - let minLocation = ref <| Point() - let max = ref 0.0 - let maxLocation = ref <| Point() - mat.MinMax(min, max, minLocation, maxLocation) - float32 !min, float32 !max - - let inline bin (x: float32) : int = - let p = int ((x - min) / (max - min) * float32 nbSamples) - if p >= nbSamples then nbSamples - 1 else p - - let data = Array.zeroCreate nbSamples - - for i in 0 .. mat.Height - 1 do - for j in 0 .. mat.Width - 1 do - let p = bin matData.[i, j] - data.[p] <- data.[p] + 1 - - { data = data; total = mat.Height * mat.Width; sum = Array.sum data; min = min; max = max } - -let histogram (values: float32 seq) (nbSamples: int) : Histogram = - let mutable min = Single.MaxValue - let mutable max = Single.MinValue - let mutable n = 0 - - for v in values do - n <- n + 1 - if v < min then min <- v - if v > max then max <- v - - let inline bin (x: float32) : int = - let p = int ((x - min) / (max - min) * float32 nbSamples) - if p >= nbSamples then nbSamples - 1 else p - - let data = Array.zeroCreate nbSamples - - for v in values do - let p = bin v - data.[p] <- data.[p] + 1 - - { data = data; total = n; sum = Array.sum data; min = min; max = max } - -let otsu (hist: Histogram) : float32 * float32 * float32 = - let mutable sumB = 0 - let mutable wB = 0 - let mutable maximum = 0.0 - let mutable level = 0 - let sum = hist.data |> Array.mapi (fun i v -> i * v |> float) |> Array.sum - - for i in 0 .. hist.data.Length - 1 do - wB <- wB + hist.data.[i] - if wB <> 0 - then - let wF = hist.total - wB - if wF <> 0 - then - sumB <- sumB + i * hist.data.[i] - let mB = (float sumB) / (float wB) - let mF = (sum - float sumB) / (float wF) - let between = (float wB) * (float wF) * (mB - mF) ** 2.; - if between >= maximum - then - level <- i - maximum <- between - - let mean1 = - let mutable sum = 0 - let mutable nb = 0 - for i in 0 .. level - 1 do - sum <- sum + i * hist.data.[i] - nb <- nb + hist.data.[i] - (sum + level * hist.data.[level] / 2) / (nb + hist.data.[level] / 2) - - let mean2 = - let mutable sum = 0 - let mutable nb = 0 - for i in level + 1 .. hist.data.Length - 1 do - sum <- sum + i * hist.data.[i] - nb <- nb + hist.data.[i] - (sum + level * hist.data.[level] / 2) / (nb + hist.data.[level] / 2) - - let toFloat l = - float32 l / float32 hist.data.Length * (hist.max - hist.min) + hist.min - - toFloat level, toFloat mean1, toFloat mean2 /// /// Remove M-adjacent pixels. It may be used after thinning. @@ -192,123 +27,6 @@ let suppressMAdjacency (img: Matrix) = then img.[i, j] <- 0uy -/// -/// Find edges of an image by using the Canny approach. -/// The thresholds are automatically defined with otsu on gradient magnitudes. -/// -/// -let findEdges (img: Image) : Matrix * Matrix * Matrix = - let w = img.Width - let h = img.Height - - use sobelKernel = - new Matrix(array2D [[ 1.0f; 0.0f; -1.0f ] - [ 2.0f; 0.0f; -2.0f ] - [ 1.0f; 0.0f; -1.0f ]]) - - let xGradient = new Matrix(img.Size) - let yGradient = new Matrix(img.Size) - CvInvoke.Filter2D(img, xGradient, sobelKernel, Point(1, 1)) - CvInvoke.Filter2D(img, yGradient, sobelKernel.Transpose(), Point(1, 1)) - - use magnitudes = new Matrix(xGradient.Size) - use angles = new Matrix(xGradient.Size) - CvInvoke.CartToPolar(xGradient, yGradient, magnitudes, angles) // Compute the magnitudes and angles. - - let thresholdHigh, thresholdLow = - let sensibilityHigh = 0.1f - let sensibilityLow = 0.0f - let threshold, _, _ = otsu (histogramMat magnitudes 300) - threshold + (sensibilityHigh * threshold), threshold - (sensibilityLow * threshold) - - // Non-maximum suppression. - use nms = new Matrix(xGradient.Size) - - let nmsData = nms.Data - let anglesData = angles.Data - let magnitudesData = magnitudes.Data - let xGradientData = xGradient.Data - let yGradientData = yGradient.Data - - for i in 0 .. h - 1 do - nmsData.[i, 0] <- 0uy - nmsData.[i, w - 1] <- 0uy - - for j in 0 .. w - 1 do - nmsData.[0, j] <- 0uy - nmsData.[h - 1, j] <- 0uy - - for i in 1 .. h - 2 do - for j in 1 .. w - 2 do - let vx = xGradientData.[i, j] - let vy = yGradientData.[i, j] - if vx <> 0.f || vy <> 0.f - then - let angle = anglesData.[i, j] - - let vx', vy' = abs vx, abs vy - let ratio2 = if vx' > vy' then vy' / vx' else vx' / vy' - let ratio1 = 1.f - ratio2 - - let mNeigbors (sign: int) : float32 = - if angle < PI / 4.f - then ratio1 * magnitudesData.[i, j + sign] + ratio2 * magnitudesData.[i + sign, j + sign] - elif angle < PI / 2.f - then ratio2 * magnitudesData.[i + sign, j + sign] + ratio1 * magnitudesData.[i + sign, j] - elif angle < 3.f * PI / 4.f - then ratio1 * magnitudesData.[i + sign, j] + ratio2 * magnitudesData.[i + sign, j - sign] - elif angle < PI - then ratio2 * magnitudesData.[i + sign, j - sign] + ratio1 * magnitudesData.[i, j - sign] - elif angle < 5.f * PI / 4.f - then ratio1 * magnitudesData.[i, j - sign] + ratio2 * magnitudesData.[i - sign, j - sign] - elif angle < 3.f * PI / 2.f - then ratio2 * magnitudesData.[i - sign, j - sign] + ratio1 * magnitudesData.[i - sign, j] - elif angle < 7.f * PI / 4.f - then ratio1 * magnitudesData.[i - sign, j] + ratio2 * magnitudesData.[i - sign, j + sign] - else ratio2 * magnitudesData.[i - sign, j + sign] + ratio1 * magnitudesData.[i, j + sign] - - let m = magnitudesData.[i, j] - if m >= thresholdLow && m > mNeigbors 1 && m > mNeigbors -1 - then - nmsData.[i, j] <- 1uy - - // suppressMConnections nms // It's not helpful for the rest of the process (ellipse detection). - - let edges = new Matrix(xGradient.Size) - let edgesData = edges.Data - - // Hysteresis thresholding. - let toVisit = Stack() - for i in 0 .. h - 1 do - for j in 0 .. w - 1 do - if nmsData.[i, j] = 1uy && magnitudesData.[i, j] >= thresholdHigh - then - nmsData.[i, j] <- 0uy - toVisit.Push(Point(j, i)) - while toVisit.Count > 0 do - let p = toVisit.Pop() - edgesData.[p.Y, p.X] <- 1uy - for i' in -1 .. 1 do - for j' in -1 .. 1 do - if i' <> 0 || j' <> 0 - then - let ni = p.Y + i' - let nj = p.X + j' - if ni >= 0 && ni < h && nj >= 0 && nj < w && nmsData.[ni, nj] = 1uy - then - nmsData.[ni, nj] <- 0uy - toVisit.Push(Point(nj, ni)) - - edges, xGradient, yGradient - -let gaussianFilter (img : Image<'TColor, 'TDepth>) (standardDeviation : float) : Image<'TColor, 'TDepth> = - let size = 2 * int (ceil (4.0 * standardDeviation)) + 1 - img.SmoothGaussian(size, size, standardDeviation, standardDeviation) - -let drawPoints (img: Image) (points: Points) (intensity: 'TDepth) = - for p in points do - img.Data.[p.Y, p.X, 0] <- intensity - type ExtremumType = | Maxima = 1 | Minima = 2 @@ -950,69 +668,3 @@ let connectedComponents (img: Image) (startPoints: List) : Po pointToCheck.Push(p) pointChecked - -let drawLine (img: Image<'TColor, 'TDepth>) (color: 'TColor) (x0: int) (y0: int) (x1: int) (y1: int) (thickness: int) = - img.Draw(LineSegment2D(Point(x0, y0), Point(x1, y1)), color, thickness); - -let drawLineF (img: Image<'TColor, 'TDepth>) (color: 'TColor) (x0: float) (y0: float) (x1: float) (y1: float) (thickness: int) = - img.Draw(LineSegment2DF(PointF(float32 x0, float32 y0), PointF(float32 x1, float32 y1)), color, thickness, CvEnum.LineType.AntiAlias); - -let drawEllipse (img: Image<'TColor, 'TDepth>) (e: Ellipse) (color: 'TColor) (alpha: float) = - if alpha >= 1.0 - then - img.Draw(Emgu.CV.Structure.Ellipse(PointF(e.Cx, e.Cy), SizeF(2.f * e.B, 2.f * e.A), e.Alpha / PI * 180.f), color, 1, CvEnum.LineType.AntiAlias) - else - let windowPosX = e.Cx - e.A - 5.f - let gapX = windowPosX - (float32 (int windowPosX)) - - let windowPosY = e.Cy - e.A - 5.f - let gapY = windowPosY - (float32 (int windowPosY)) - - let roi = Rectangle(int windowPosX, int windowPosY, 2.f * (e.A + 5.f) |> int, 2.f * (e.A + 5.f) |> int) - - img.ROI <- roi - if roi = img.ROI // We do not display ellipses touching the edges (FIXME) - then - use i = new Image<'TColor, 'TDepth>(img.ROI.Size) - i.Draw(Emgu.CV.Structure.Ellipse(PointF(e.A + 5.f + gapX, e.A + 5.f + gapY), SizeF(2.f * e.B, 2.f * e.A), e.Alpha / PI * 180.f), color, 1, CvEnum.LineType.AntiAlias) - CvInvoke.AddWeighted(img, 1.0, i, alpha, 0.0, img) - img.ROI <- Rectangle.Empty - -let drawEllipses (img: Image<'TColor, 'TDepth>) (ellipses: Ellipse list) (color: 'TColor) (alpha: float) = - List.iter (fun e -> drawEllipse img e color alpha) ellipses - -let rngCell = System.Random() -let drawCell (img: Image) (drawCellContent: bool) (c: Cell) = - if drawCellContent - then - let colorB = rngCell.Next(20, 70) - let colorG = rngCell.Next(20, 70) - let colorR = rngCell.Next(20, 70) - - for y in 0 .. c.elements.Height - 1 do - for x in 0 .. c.elements.Width - 1 do - if c.elements.[y, x] > 0uy - then - let dx, dy = c.center.X - c.elements.Width / 2, c.center.Y - c.elements.Height / 2 - let b = img.Data.[y + dy, x + dx, 0] |> int - let g = img.Data.[y + dy, x + dx, 1] |> int - let r = img.Data.[y + dy, x + dx, 2] |> int - img.Data.[y + dy, x + dx, 0] <- if b + colorB > 255 then 255uy else byte (b + colorB) - img.Data.[y + dy, x + dx, 1] <- if g + colorG > 255 then 255uy else byte (g + colorG) - img.Data.[y + dy, x + dx, 2] <- if r + colorR > 255 then 255uy else byte (r + colorR) - - let crossColor, crossColor2 = - match c.cellClass with - | HealthyRBC -> Bgr(255., 0., 0.), Bgr(255., 255., 255.) - | InfectedRBC -> Bgr(0., 0., 255.), Bgr(120., 120., 255.) - | Peculiar -> Bgr(0., 0., 0.), Bgr(80., 80., 80.) - - drawLine img crossColor2 (c.center.X - 3) c.center.Y (c.center.X + 3) c.center.Y 2 - drawLine img crossColor2 c.center.X (c.center.Y - 3) c.center.X (c.center.Y + 3) 2 - - drawLine img crossColor (c.center.X - 3) c.center.Y (c.center.X + 3) c.center.Y 1 - drawLine img crossColor c.center.X (c.center.Y - 3) c.center.X (c.center.Y + 3) 1 - - -let drawCells (img: Image) (drawCellContent: bool) (cells: Cell list) = - List.iter (fun c -> drawCell img drawCellContent c) cells \ No newline at end of file diff --git a/Parasitemia/ParasitemiaCore/ImgTools/Otsu.fs b/Parasitemia/ParasitemiaCore/ImgTools/Otsu.fs new file mode 100644 index 0000000..6b8ee53 --- /dev/null +++ b/Parasitemia/ParasitemiaCore/ImgTools/Otsu.fs @@ -0,0 +1,47 @@ +module ParasitemiaCore.Otsu + +open Histogram + +let otsu (hist: Histogram) : float32 * float32 * float32 = + let mutable sumB = 0 + let mutable wB = 0 + let mutable maximum = 0.0 + let mutable level = 0 + let sum = hist.data |> Array.mapi (fun i v -> i * v |> float) |> Array.sum + + for i in 0 .. hist.data.Length - 1 do + wB <- wB + hist.data.[i] + if wB <> 0 + then + let wF = hist.total - wB + if wF <> 0 + then + sumB <- sumB + i * hist.data.[i] + let mB = (float sumB) / (float wB) + let mF = (sum - float sumB) / (float wF) + let between = (float wB) * (float wF) * (mB - mF) ** 2.; + if between >= maximum + then + level <- i + maximum <- between + + let mean1 = + let mutable sum = 0 + let mutable nb = 0 + for i in 0 .. level - 1 do + sum <- sum + i * hist.data.[i] + nb <- nb + hist.data.[i] + (sum + level * hist.data.[level] / 2) / (nb + hist.data.[level] / 2) + + let mean2 = + let mutable sum = 0 + let mutable nb = 0 + for i in level + 1 .. hist.data.Length - 1 do + sum <- sum + i * hist.data.[i] + nb <- nb + hist.data.[i] + (sum + level * hist.data.[level] / 2) / (nb + hist.data.[level] / 2) + + let toFloat l = + float32 l / float32 hist.data.Length * (hist.max - hist.min) + hist.min + + toFloat level, toFloat mean1, toFloat mean2 \ No newline at end of file diff --git a/Parasitemia/ParasitemiaCore/MainAnalysis.fs b/Parasitemia/ParasitemiaCore/MainAnalysis.fs index 9c85c7b..1d2ced3 100644 --- a/Parasitemia/ParasitemiaCore/MainAnalysis.fs +++ b/Parasitemia/ParasitemiaCore/MainAnalysis.fs @@ -12,6 +12,7 @@ open Emgu.CV.Structure open Logger open Utils +open Morpho open ImgTools open Config open Types @@ -88,7 +89,7 @@ let doAnalysis (img: Image) (name: string) (config: Config) (reportPr logTimeWithName "Parasites segmentation" (fun () -> reportWithVal 40 (ParasitesMarker.find img_parasites_filtered config)) let! edges, xGradient, yGradient = logTimeWithName "Finding edges" (fun () -> - let edges, xGradient, yGradient = findEdges img_RBC_filtered + let edges, xGradient, yGradient = Edges.find img_RBC_filtered removeArea edges (config.RBCRadius.Pixel ** 2.f / 50.f |> int) reportWithVal 50 (edges, xGradient, yGradient)) @@ -109,43 +110,43 @@ let doAnalysis (img: Image) (name: string) (config: Config) (reportPr let buildFileName postfix = System.IO.Path.Combine(dirPath, name + postfix) - saveMat (edges * 255.0) (buildFileName " - edges.png") + IO.saveMat (edges * 255.0) (buildFileName " - edges.png") - saveImg parasites.darkStain (buildFileName " - parasites - dark stain.png") - saveImg parasites.parasite (buildFileName " - parasites - stain.png") - saveImg parasites.nucleus (buildFileName " - parasites - infection.png") + IO.saveImg parasites.darkStain (buildFileName " - parasites - dark stain.png") + IO.saveImg parasites.parasite (buildFileName " - parasites - stain.png") + IO.saveImg parasites.nucleus (buildFileName " - parasites - infection.png") let imgAllEllipses = img.Copy() - drawEllipses imgAllEllipses matchingEllipses.Ellipses (Bgr(255.0, 255.0, 255.0)) 0.04 - saveImg imgAllEllipses (buildFileName " - ellipses - all.png") + Drawing.drawEllipses imgAllEllipses matchingEllipses.Ellipses (Bgr(255.0, 255.0, 255.0)) 0.04 + IO.saveImg imgAllEllipses (buildFileName " - ellipses - all.png") let imgEllipses = img_RBC_filtered.Convert() - drawEllipses imgEllipses prunedEllipses (Bgr(0.0, 240.0, 240.0)) 1.0 - saveImg imgEllipses (buildFileName " - ellipses.png") + Drawing.drawEllipses imgEllipses prunedEllipses (Bgr(0.0, 240.0, 240.0)) 1.0 + IO.saveImg imgEllipses (buildFileName " - ellipses.png") let imgCells = img.Copy() - drawCells imgCells false cells - saveImg imgCells (buildFileName " - cells.png") + Drawing.drawCells imgCells false cells + IO.saveImg imgCells (buildFileName " - cells.png") let imgCells' = img.Copy() - drawCells imgCells' true cells - saveImg imgCells' (buildFileName " - cells - full.png") + Drawing.drawCells imgCells' true cells + IO.saveImg imgCells' (buildFileName " - cells - full.png") let filteredGreenMaxima = gaussianFilter img_RBC config.LPFStandardDeviationRBC for m in findMaxima filteredGreenMaxima do - drawPoints filteredGreenMaxima m 255.f - saveImg filteredGreenMaxima (buildFileName " - filtered - maxima.png") + Drawing.drawPoints filteredGreenMaxima m 255.f + IO.saveImg filteredGreenMaxima (buildFileName " - filtered - maxima.png") - saveImg img_RBC_filtered (buildFileName " - filtered.png") - saveImg imgWhitoutParasite (buildFileName " - filtered closed stain.png") - saveImg imgWithoutNucleus (buildFileName " - filtered closed infection.png") + IO.saveImg img_RBC_filtered (buildFileName " - filtered.png") + IO.saveImg imgWhitoutParasite (buildFileName " - filtered closed stain.png") + IO.saveImg imgWithoutNucleus (buildFileName " - filtered closed infection.png") - saveImg img_RBC (buildFileName " - source - RBC.png") - saveImg img_parasites (buildFileName " - source - parasites.png") + IO.saveImg img_RBC (buildFileName " - source - RBC.png") + IO.saveImg img_parasites (buildFileName " - source - parasites.png") - saveImg (normalize img_float.[2] 255.) (buildFileName " - source - red.png") - saveImg (normalize img_float.[1] 255.) (buildFileName " - source - green.png") - saveImg (normalize img_float.[0] 255.) (buildFileName " - source - blue.png") + IO.saveImg (normalize img_float.[2] 255.) (buildFileName " - source - red.png") + IO.saveImg (normalize img_float.[1] 255.) (buildFileName " - source - green.png") + IO.saveImg (normalize img_float.[0] 255.) (buildFileName " - source - blue.png") | _ -> () return cells } diff --git a/Parasitemia/ParasitemiaCore/ParasitemiaCore.fsproj b/Parasitemia/ParasitemiaCore/ParasitemiaCore.fsproj index 58886c0..6a9eb75 100644 --- a/Parasitemia/ParasitemiaCore/ParasitemiaCore.fsproj +++ b/Parasitemia/ParasitemiaCore/ParasitemiaCore.fsproj @@ -51,7 +51,13 @@ - + + + + + + + diff --git a/Parasitemia/ParasitemiaCore/ParasitesMarker.fs b/Parasitemia/ParasitemiaCore/ParasitesMarker.fs index 6c25d04..0d7c7ae 100644 --- a/Parasitemia/ParasitemiaCore/ParasitesMarker.fs +++ b/Parasitemia/ParasitemiaCore/ParasitesMarker.fs @@ -7,6 +7,9 @@ open Emgu.CV open Emgu.CV.Structure open Utils +open Histogram +open Otsu +open Morpho open ImgTools type Result = { diff --git a/Parasitemia/ParasitemiaCore/UnitsOfMeasure.fs b/Parasitemia/ParasitemiaCore/UnitsOfMeasure.fs index 0d84d16..36c9495 100644 --- a/Parasitemia/ParasitemiaCore/UnitsOfMeasure.fs +++ b/Parasitemia/ParasitemiaCore/UnitsOfMeasure.fs @@ -1,9 +1,9 @@ module ParasitemiaCore.UnitsOfMeasure -[] type px +[] type px // Pixel. [] type μm [] type inch -[] type ppi = px / inch +[] type ppi = px / inch // Pixel per inch. let μmInchRatio = 25.4e3<μm/inch> diff --git a/Parasitemia/ParasitemiaUI/XAML/AboutWindow.xaml b/Parasitemia/ParasitemiaUI/XAML/AboutWindow.xaml index 3cb0de4..e35d234 100644 --- a/Parasitemia/ParasitemiaUI/XAML/AboutWindow.xaml +++ b/Parasitemia/ParasitemiaUI/XAML/AboutWindow.xaml @@ -3,7 +3,7 @@ xmlns:d="http://schemas.microsoft.com/expression/blend/2008" xmlns:mc="http://schemas.openxmlformats.org/markup-compatibility/2006" mc:Ignorable="d" - x:Name="AboutWindow" Height="200.969" Width="282.313" MinHeight="100" MinWidth="100" Title="About" Icon="pack://application:,,,/Resources/icon.ico"> + x:Name="AboutWindow" Height="220" Width="280" MinHeight="220" MinWidth="280" Title="About" Icon="pack://application:,,,/Resources/icon.ico" ResizeMode="NoResize"> -- 2.45.2