X-Git-Url: http://git.euphorik.ch/?p=master-thesis.git;a=blobdiff_plain;f=Parasitemia%2FParasitemia%2FEEOver.fs;h=528632ad29fc856a6fb78d03015611dbc0e9b3a0;hp=ef21051ddbfbb5159e9d1d9c6f2af9447610a9bd;hb=53440e757b4a4ab2a81b0f6a5dd1a2002c0133ba;hpb=962440dfe579569ba7e3ebba6f878cb575d1fc02 diff --git a/Parasitemia/Parasitemia/EEOver.fs b/Parasitemia/Parasitemia/EEOver.fs index ef21051..528632a 100644 --- a/Parasitemia/Parasitemia/EEOver.fs +++ b/Parasitemia/Parasitemia/EEOver.fs @@ -437,33 +437,33 @@ let private biquadroots (p: float[]) (r: float[,]) = p.[k] <- p.[k] / p.[0] p.[0] <- 1.0 let e = 0.25 * p.[1] - let b = ref (2.0 * e) - let c = ref (!b ** 2.0) - let mutable d = 0.75 * !c - b := p.[3] + !b *(!c - p.[2]) + let mutable b = 2.0 * e + let mutable c = b ** 2.0 + let mutable d = 0.75 * c + b <- p.[3] + b *(c - p.[2]) let mutable a = p.[2] - d - c := p.[4] + e * (e * a - p.[3]) + c <- p.[4] + e * (e * a - p.[3]) a <- a - d - let quadExecuted = ref false + let mutable quadExecuted = false let quad () = - if not !quadExecuted + if not quadExecuted then - p.[2] <- !c / !b + p.[2] <- c / b quadroots p r for k in 1..2 do for j in 1..2 do r.[j, k+2] <- r.[j, k] p.[1] <- -p.[1] - p.[2] <- !b + p.[2] <- b quadroots p r for k in 1..4 do r.[1,k] <- r.[1,k] - e - quadExecuted := true + quadExecuted <- true p.[1] <- 0.5 * a - p.[2] <- (p.[1] * p.[1] - !c) * 0.25 - p.[3] <- !b * !b / -64.0 + p.[2] <- (p.[1] * p.[1] - c) * 0.25 + p.[3] <- b * b / -64.0 if p.[3] < 0.0 then cubicroots p r @@ -473,43 +473,43 @@ let private biquadroots (p: float[]) (r: float[,]) = then d <- r.[1, k] * 4.0 a <- a + d - if a >= 0.0 && !b >= 0.0 + if a >= 0.0 && b >= 0.0 then p.[1] <- sqrt d - elif a <= 0.0 && !b <= 0.0 + elif a <= 0.0 && b <= 0.0 then p.[1] <- sqrt d else p.[1] <- -(sqrt d) - b := 0.5 * (a + !b / p.[1]) + b <- 0.5 * (a + b / p.[1]) quad () k <- 4 k <- k + 1 - if not !quadExecuted && p.[2] < 0.0 + if not quadExecuted && p.[2] < 0.0 then - b := sqrt !c - d <- !b + !b - a + b <- sqrt c + d <- b + b - a p.[1] <- 0.0 if d > 0.0 then p.[1] <- sqrt d - elif not !quadExecuted + elif not quadExecuted then if p.[1] > 0.0 then - b := (sqrt p.[2]) * 2.0 + p.[1] + b <- (sqrt p.[2]) * 2.0 + p.[1] else - b := -(sqrt p.[2]) * 2.0 + p.[1] + b <- -(sqrt p.[2]) * 2.0 + p.[1] - if !b <> 0.0 + if b <> 0.0 then p.[1] <- 0.0 else for k in 1..4 do r.[1, k] <- -e r.[2, k] <- 0.0 - quadExecuted := true + quadExecuted <- true quad () @@ -610,15 +610,16 @@ let EEOverlapArea (e1: Types.Ellipse) (e2: Types.Ellipse) : (float * float[] * f printf "nroots = %d\n" nroots #endif - let ychk = [| - for i in 1 .. nroots do - if abs r.[2, i] < EPS - then - yield r.[1, i] * b1 + let ychk = Array.init nroots (fun _ -> Double.MaxValue) + let mutable nychk = 0 + for i in 1 .. nroots do + if abs r.[2, i] < EPS + then + ychk.[nychk] <- r.[1, i] * b1 + nychk <- nychk + 1 #if DEBUG_LOG - printf "ROOT is Real, i=%d --> %f (B1=%f)\n" i r.[1, i] b1 + printf "ROOT is Real, i=%d --> %f (B1=%f)\n" i r.[1, i] b1 #endif - |] Array.sortInPlace ychk #if DEBUG_LOG @@ -627,7 +628,6 @@ let EEOverlapArea (e1: Types.Ellipse) (e2: Types.Ellipse) : (float * float[] * f printf "\t j=%d, ychk=%f\n" j ychk.[j] #endif - let nychk = Array.length ychk let mutable nintpts = 0 let xint = Array.zeroCreate 4