X-Git-Url: http://git.euphorik.ch/?p=master-thesis.git;a=blobdiff_plain;f=Parasitemia%2FParasitemia%2FImgTools.fs;h=e9cbbe9ee0f59f7a4f54c557ff345ca2f20be6df;hp=f7e5a4af45a969c27d5fe69abdcc5240ffb86d5b;hb=b070295cf67b2025164a34b6594e84f0d771cdc9;hpb=5b68d946369f998865d2dd330fd3b374b2b9a0ad diff --git a/Parasitemia/Parasitemia/ImgTools.fs b/Parasitemia/Parasitemia/ImgTools.fs index f7e5a4a..e9cbbe9 100644 --- a/Parasitemia/Parasitemia/ImgTools.fs +++ b/Parasitemia/Parasitemia/ImgTools.fs @@ -3,11 +3,13 @@ open System open System.Drawing open System.Collections.Generic +open System.Linq open Emgu.CV open Emgu.CV.Structure open Utils +open Heap // Normalize image values between 0uy and 255uy. let normalizeAndConvert (img: Image) : Image = @@ -18,10 +20,266 @@ let normalizeAndConvert (img: Image) : Image = img.MinMax(min, max, minLocation, maxLocation) ((img - (!min).[0]) / ((!max).[0] - (!min).[0]) * 255.0).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) + +let findMaxima (img: Image) : IEnumerable> = + use suppress = new Image(img.Size) + let w = img.Width + let h = img.Height + + let imgData = img.Data + let suppressData = suppress.Data + + let result = List>() + + let flood (start: Point) : List> = + let sameLevelToCheck = Stack() + let betterLevelToCheck = Stack() + betterLevelToCheck.Push(start) + + let result' = List>() + + while betterLevelToCheck.Count > 0 do + let p = betterLevelToCheck.Pop() + if suppressData.[p.Y, p.X, 0] = 0uy + then + suppressData.[p.Y, p.X, 0] <- 1uy + sameLevelToCheck.Push(p) + let current = List() + + let mutable betterExists = false + + while sameLevelToCheck.Count > 0 do + let p' = sameLevelToCheck.Pop() + let currentLevel = imgData.[p'.Y, p'.X, 0] + current.Add(p') |> ignore + for i in -1 .. 1 do + for j in -1 .. 1 do + if i <> 0 || j <> 0 + then + let ni = i + p'.Y + let nj = j + p'.X + if ni >= 0 && ni < h && nj >= 0 && nj < w + then + let level = imgData.[ni, nj, 0] + let notSuppressed = suppressData.[ni, nj, 0] = 0uy + + if level = currentLevel && notSuppressed + then + suppressData.[ni, nj, 0] <- 1uy + sameLevelToCheck.Push(Point(nj, ni)) + elif level > currentLevel + then + betterExists <- true + if notSuppressed + then + betterLevelToCheck.Push(Point(nj, ni)) + + if not betterExists + then + result'.Add(current) + result' + + for i in 0 .. h - 1 do + for j in 0 .. w - 1 do + let maxima = flood (Point(j, i)) + if maxima.Count > 0 + then result.AddRange(maxima) + + result.Select(fun l -> HashSet(l)) + + +type PriorityQueue () = + let q = List>() // TODO: Check performance with an HasSet + let mutable highest = -1 // Value of the first elements of 'q'. + + member this.Next () : byte * Point = + if this.IsEmpty + then + invalidOp "Queue is empty" + else + let l = q.[0] + let next = l.First() + l.Remove(next) |> ignore + let value = byte highest + if l.Count = 0 + then + q.RemoveAt(0) + highest <- highest - 1 + while q.Count > 0 && q.[0] = null do + q.RemoveAt(0) + highest <- highest - 1 + value, next + + member this.Max = + highest |> byte + + (*member this.UnionWith (other: PriorityQueue) = + while not other.IsEmpty do + let p, v = other.Next + this.Add p v*) + + member this.Add (value: byte) (p: Point) = + let vi = int value + + if this.IsEmpty + then + highest <- int value + q.Insert(0, null) + elif vi > highest + then + for i in highest .. vi - 1 do + q.Insert(0, null) + highest <- vi + elif highest - vi >= q.Count + then + for i in 0 .. highest - vi - q.Count do + q.Add(null) + + let pos = highest - vi + if q.[pos] = null + then + q.[pos] <- HashSet([p]) + else + q.[pos].Add(p) |> ignore + + + member this.IsEmpty = + q.Count = 0 + + member this.Clear () = + while highest >= 0 do + q.[highest] <- null + highest <- highest - 1 + + + +type MaximaState = Uncertain | Validated | TooBig +type Maxima = { + elements : HashSet + mutable intensity: byte option + mutable state: MaximaState } + + +let areaOpen (img: Image) (area: int) = + let w = img.Width + let h = img.Height + + let maxima = findMaxima img |> Seq.map (fun m -> { elements = m; intensity = None; state = Uncertain }) |> List.ofSeq + let toValidated = Stack(maxima) + + while toValidated.Count > 0 do + let m = toValidated.Pop() + if m.elements.Count <= area + then + let queue = + let q = PriorityQueue() + let firstElements = HashSet() + for p in m.elements do + for i in -1 .. 1 do + for j in -1 .. 1 do + if i <> 0 || j <> 0 + then + let ni = i + p.Y + let nj = j + p.X + let p' = Point(nj, ni) + if ni >= 0 && ni < h && nj >= 0 && nj < w && not (m.elements.Contains(p')) && not (firstElements.Contains(p')) + then + firstElements.Add(p') |> ignore + q.Add (img.Data.[ni, nj, 0]) p' + q + + let mutable intensity = queue.Max + let nextElements = HashSet() + + let mutable stop = false + while not stop do + let intensity', p = queue.Next () + + if intensity' = intensity // The intensity doesn't change. + then + if m.elements.Count + nextElements.Count + 1 > area + then + m.state <- Validated + m.intensity <- Some intensity + stop <- true + else + nextElements.Add(p) |> ignore + elif intensity' < intensity + then + m.elements.UnionWith(nextElements) + if m.elements.Count = area + then + m.state <- Validated + m.intensity <- Some (intensity') + stop <- true + else + intensity <- intensity' + nextElements.Clear() + nextElements.Add(p) |> ignore + else // i' > i + seq { + for m' in maxima do + if m' <> m && m'.elements.Contains(p) then + if m'.elements.Count + m.elements.Count <= area + then + m'.state <- Uncertain + m'.elements.UnionWith(m.elements) + if not <| toValidated.Contains m' // FIXME: Maybe use state instead of scanning the whole list. + then + toValidated.Push(m') + stop <- true + yield false + } |> Seq.forall id |> ignore + + if not stop + then + m.state <- Validated + m.intensity <- Some (intensity) + stop <- true + + if not stop + then + for i in -1 .. 1 do + for j in -1 .. 1 do + if i <> 0 || j <> 0 + then + let ni = i + p.Y + let nj = j + p.X + let p' = Point(nj, ni) + if ni < 0 || ni >= h || nj < 0 || nj >= w + then + m.state <- Validated + m.intensity <- Some (intensity) + stop <- true + elif not (m.elements.Contains(p')) && not (nextElements.Contains(p')) + then + queue.Add (img.Data.[ni, nj, 0]) p' + + if queue.IsEmpty + then + if m.elements.Count + nextElements.Count <= area + then + m.state <- Validated + m.intensity <- Some intensity' + m.elements.UnionWith(nextElements) + stop <- true + + for m in maxima do + if m.state = Validated + then + match m.intensity with + | Some i -> + for p in m.elements do + img.Data.[p.Y, p.X, 0] <- i + | _ -> () + () + + // Zhang and Suen algorithm. // Modify 'mat' in place. let thin (mat: Matrix) = @@ -73,6 +331,7 @@ let thin (mat: Matrix) = data2 <- tmp +// FIXME: replace by a queue or stack. let pop (l: List<'a>) : 'a = let n = l.[l.Count - 1] l.RemoveAt(l.Count - 1) @@ -156,41 +415,37 @@ let saveMat (mat: Matrix<'TDepth>) (filepath: string) = mat.CopyTo(img) saveImg img filepath -let drawLine (img: Image<'TColor, 'TDepth>) (color: 'TColor) (x0: int) (y0: int) (x1: int) (y1: int) = - img.Draw(LineSegment2D(Point(x0, y0), Point(x1, y1)), color, 1); +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) = - let x0, y0, x1, y1 = roundInt(x0), roundInt(y0), roundInt(x1), roundInt(y1) - drawLine img color x0 y0 x1 y1 +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: Types.Ellipse) (color: 'TColor) = - let cosAlpha = cos e.Alpha - let sinAlpha = sin e.Alpha +let drawEllipse (img: Image<'TColor, 'TDepth>) (e: Types.Ellipse) (color: 'TColor) (alpha: float) = - let mutable x0 = 0.0 - let mutable y0 = 0.0 - let mutable first_iteration = true + if alpha >= 1.0 + then + img.Draw(Ellipse(PointF(float32 e.Cx, float32 e.Cy), SizeF(2. * e.B |> float32, 2. * e.A |> float32), float32 <| e.Alpha / Math.PI * 180.), color, 1, CvEnum.LineType.AntiAlias) + else + let windowPosX = e.Cx - e.A - 5.0 + let gapX = windowPosX - (float (int windowPosX)) - let n = 40 - let thetaIncrement = 2.0 * Math.PI / (float n) + let windowPosY = e.Cy - e.A - 5.0 + let gapY = windowPosY - (float (int windowPosY)) - for theta in 0.0 .. thetaIncrement .. 2.0 * Math.PI do - let cosTheta = cos theta - let sinTheta = sin theta - let x = e.Cx + cosAlpha * e.A * cosTheta - sinAlpha * e.B * sinTheta - let y = e.Cy + sinAlpha * e.A * cosTheta + cosAlpha * e.B * sinTheta + let roi = Rectangle(int windowPosX, int windowPosY, 2. * (e.A + 5.0) |> int, 2.* (e.A + 5.0) |> int) - if not first_iteration + img.ROI <- roi + if roi = img.ROI // We do not display ellipses touching the edges (FIXME) then - drawLineF img color x0 y0 x y - else - first_iteration <- false + use i = new Image<'TColor, 'TDepth>(img.ROI.Size) + i.Draw(Ellipse(PointF(float32 <| (e.A + 5. + gapX) , float32 <| (e.A + 5. + gapY)), SizeF(2. * e.B |> float32, 2. * e.A |> float32), float32 <| e.Alpha / Math.PI * 180.), color, 1, CvEnum.LineType.AntiAlias) + CvInvoke.AddWeighted(img, 1.0, i, alpha, 0.0, img) + img.ROI <- Rectangle.Empty - x0 <- x - y0 <- y -let drawEllipses (img: Image<'TColor, 'TDepth>) (ellipses: Types.Ellipse list) (color: 'TColor) = - List.iter (fun e -> drawEllipse img e color) ellipses +let drawEllipses (img: Image<'TColor, 'TDepth>) (ellipses: Types.Ellipse list) (color: 'TColor) (alpha: float) = + List.iter (fun e -> drawEllipse img e color alpha) ellipses let rngCell = System.Random() @@ -213,13 +468,18 @@ let drawCell (img: Image) (drawCellContent: bool) (c: Types.Cell) = 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 = match c.cellClass with - | Types.HealthyRBC -> Bgr(255.0, 0.0, 0.0) - | Types.InfectedRBC -> Bgr(0.0, 0.0, 255.0) - | Types.Peculiar -> Bgr(0.0, 0.0, 0.0) + let crossColor, crossColor2 = + match c.cellClass with + | Types.HealthyRBC -> Bgr(255., 0., 0.), Bgr(255., 255., 255.) + | Types.InfectedRBC -> Bgr(0., 0., 255.), Bgr(120., 120., 255.) + | Types.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 - drawLine img crossColor (c.center.X - 3) c.center.Y (c.center.X + 3) c.center.Y - drawLine img crossColor c.center.X (c.center.Y - 3) c.center.X (c.center.Y + 3) let drawCells (img: Image) (drawCellContent: bool) (cells: Types.Cell list) = List.iter (fun c -> drawCell img drawCellContent c) cells \ No newline at end of file