1// Copyright 2016 The Go Authors. All rights reserved.
2// Use of this source code is governed by a BSD-style
3// license that can be found in the LICENSE file.
4
5package vector
6
7// This file contains a floating point math implementation of the vector
8// graphics rasterizer.
9
10import (
11 "math"
12)
13
14func floatingMax(x, y float32) float32 {
15 if x > y {
16 return x
17 }
18 return y
19}
20
21func floatingMin(x, y float32) float32 {
22 if x < y {
23 return x
24 }
25 return y
26}
27
28func floatingFloor(x float32) int32 { return int32(math.Floor(float64(x))) }
29func floatingCeil(x float32) int32 { return int32(math.Ceil(float64(x))) }
30
31func (z *Rasterizer) floatingLineTo(bx, by float32) {
32 ax, ay := z.penX, z.penY
33 z.penX, z.penY = bx, by
34 dir := float32(1)
35 if ay > by {
36 dir, ax, ay, bx, by = -1, bx, by, ax, ay
37 }
38 // Horizontal line segments yield no change in coverage. Almost horizontal
39 // segments would yield some change, in ideal math, but the computation
40 // further below, involving 1 / (by - ay), is unstable in floating point
41 // math, so we treat the segment as if it was perfectly horizontal.
42 if by-ay <= 0.000001 {
43 return
44 }
45 dxdy := (bx - ax) / (by - ay)
46
47 x := ax
48 y := floatingFloor(ay)
49 yMax := floatingCeil(by)
50 if yMax > int32(z.size.Y) {
51 yMax = int32(z.size.Y)
52 }
53 width := int32(z.size.X)
54
55 for ; y < yMax; y++ {
56 dy := floatingMin(float32(y+1), by) - floatingMax(float32(y), ay)
57
58 // The "float32" in expressions like "float32(foo*bar)" here and below
59 // look redundant, since foo and bar already have type float32, but are
60 // explicit in order to disable the compiler's Fused Multiply Add (FMA)
61 // instruction selection, which can improve performance but can result
62 // in different rounding errors in floating point computations.
63 //
64 // This package aims to have bit-exact identical results across all
65 // GOARCHes, and across pure Go code and assembly, so it disables FMA.
66 //
67 // See the discussion at
68 // https://groups.google.com/d/topic/golang-dev/Sti0bl2xUXQ/discussion
69 xNext := x + float32(dy*dxdy)
70 if y < 0 {
71 x = xNext
72 continue
73 }
74 buf := z.bufF32[y*width:]
75 d := float32(dy * dir)
76 x0, x1 := x, xNext
77 if x > xNext {
78 x0, x1 = x1, x0
79 }
80 x0i := floatingFloor(x0)
81 x0Floor := float32(x0i)
82 x1i := floatingCeil(x1)
83 x1Ceil := float32(x1i)
84
85 if x1i <= x0i+1 {
86 xmf := float32(0.5*(x+xNext)) - x0Floor
87 if i := clamp(x0i+0, width); i < uint(len(buf)) {
88 buf[i] += d - float32(d*xmf)
89 }
90 if i := clamp(x0i+1, width); i < uint(len(buf)) {
91 buf[i] += float32(d * xmf)
92 }
93 } else {
94 s := 1 / (x1 - x0)
95 x0f := x0 - x0Floor
96 oneMinusX0f := 1 - x0f
97 a0 := float32(0.5 * s * oneMinusX0f * oneMinusX0f)
98 x1f := x1 - x1Ceil + 1
99 am := float32(0.5 * s * x1f * x1f)
100
101 if i := clamp(x0i, width); i < uint(len(buf)) {
102 buf[i] += float32(d * a0)
103 }
104
105 if x1i == x0i+2 {
106 if i := clamp(x0i+1, width); i < uint(len(buf)) {
107 buf[i] += float32(d * (1 - a0 - am))
108 }
109 } else {
110 a1 := float32(s * (1.5 - x0f))
111 if i := clamp(x0i+1, width); i < uint(len(buf)) {
112 buf[i] += float32(d * (a1 - a0))
113 }
114 dTimesS := float32(d * s)
115 for xi := x0i + 2; xi < x1i-1; xi++ {
116 if i := clamp(xi, width); i < uint(len(buf)) {
117 buf[i] += dTimesS
118 }
119 }
120 a2 := a1 + float32(s*float32(x1i-x0i-3))
121 if i := clamp(x1i-1, width); i < uint(len(buf)) {
122 buf[i] += float32(d * (1 - a2 - am))
123 }
124 }
125
126 if i := clamp(x1i, width); i < uint(len(buf)) {
127 buf[i] += float32(d * am)
128 }
129 }
130
131 x = xNext
132 }
133}
134
135const (
136 // almost256 scales a floating point value in the range [0, 1] to a uint8
137 // value in the range [0x00, 0xff].
138 //
139 // 255 is too small. Floating point math accumulates rounding errors, so a
140 // fully covered src value that would in ideal math be float32(1) might be
141 // float32(1-ε), and uint8(255 * (1-ε)) would be 0xfe instead of 0xff. The
142 // uint8 conversion rounds to zero, not to nearest.
143 //
144 // 256 is too big. If we multiplied by 256, below, then a fully covered src
145 // value of float32(1) would translate to uint8(256 * 1), which can be 0x00
146 // instead of the maximal value 0xff.
147 //
148 // math.Float32bits(almost256) is 0x437fffff.
149 almost256 = 255.99998
150
151 // almost65536 scales a floating point value in the range [0, 1] to a
152 // uint16 value in the range [0x0000, 0xffff].
153 //
154 // math.Float32bits(almost65536) is 0x477fffff.
155 almost65536 = almost256 * 256
156)
157
158func floatingAccumulateOpOver(dst []uint8, src []float32) {
159 // Sanity check that len(dst) >= len(src).
160 if len(dst) < len(src) {
161 return
162 }
163
164 acc := float32(0)
165 for i, v := range src {
166 acc += v
167 a := acc
168 if a < 0 {
169 a = -a
170 }
171 if a > 1 {
172 a = 1
173 }
174 // This algorithm comes from the standard library's image/draw package.
175 dstA := uint32(dst[i]) * 0x101
176 maskA := uint32(almost65536 * a)
177 outA := dstA*(0xffff-maskA)/0xffff + maskA
178 dst[i] = uint8(outA >> 8)
179 }
180}
181
182func floatingAccumulateOpSrc(dst []uint8, src []float32) {
183 // Sanity check that len(dst) >= len(src).
184 if len(dst) < len(src) {
185 return
186 }
187
188 acc := float32(0)
189 for i, v := range src {
190 acc += v
191 a := acc
192 if a < 0 {
193 a = -a
194 }
195 if a > 1 {
196 a = 1
197 }
198 dst[i] = uint8(almost256 * a)
199 }
200}
201
202func floatingAccumulateMask(dst []uint32, src []float32) {
203 // Sanity check that len(dst) >= len(src).
204 if len(dst) < len(src) {
205 return
206 }
207
208 acc := float32(0)
209 for i, v := range src {
210 acc += v
211 a := acc
212 if a < 0 {
213 a = -a
214 }
215 if a > 1 {
216 a = 1
217 }
218 dst[i] = uint32(almost65536 * a)
219 }
220}