convolve.go 6.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274
  1. // Copyright 2011 The Graphics-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. package convolve
  5. import (
  6. "errors"
  7. "fmt"
  8. "image"
  9. "image/draw"
  10. "math"
  11. )
  12. // clamp clamps x to the range [x0, x1].
  13. func clamp(x, x0, x1 float64) float64 {
  14. if x < x0 {
  15. return x0
  16. }
  17. if x > x1 {
  18. return x1
  19. }
  20. return x
  21. }
  22. // Kernel is a square matrix that defines a convolution.
  23. type Kernel interface {
  24. // Weights returns the square matrix of weights in row major order.
  25. Weights() []float64
  26. }
  27. // SeparableKernel is a linearly separable, square convolution kernel.
  28. // X and Y are the per-axis weights. Each slice must be the same length, and
  29. // have an odd length. The middle element of each slice is the weight for the
  30. // central pixel. For example, the horizontal Sobel kernel is:
  31. // sobelX := &SeparableKernel{
  32. // X: []float64{-1, 0, +1},
  33. // Y: []float64{1, 2, 1},
  34. // }
  35. type SeparableKernel struct {
  36. X, Y []float64
  37. }
  38. func (k *SeparableKernel) Weights() []float64 {
  39. n := len(k.X)
  40. w := make([]float64, n*n)
  41. for y := 0; y < n; y++ {
  42. for x := 0; x < n; x++ {
  43. w[y*n+x] = k.X[x] * k.Y[y]
  44. }
  45. }
  46. return w
  47. }
  48. // fullKernel is a square convolution kernel.
  49. type fullKernel []float64
  50. func (k fullKernel) Weights() []float64 { return k }
  51. func kernelSize(w []float64) (size int, err error) {
  52. size = int(math.Sqrt(float64(len(w))))
  53. if size*size != len(w) {
  54. return 0, errors.New("graphics: kernel is not square")
  55. }
  56. if size%2 != 1 {
  57. return 0, errors.New("graphics: kernel size is not odd")
  58. }
  59. return size, nil
  60. }
  61. // NewKernel returns a square convolution kernel.
  62. func NewKernel(w []float64) (Kernel, error) {
  63. if _, err := kernelSize(w); err != nil {
  64. return nil, err
  65. }
  66. return fullKernel(w), nil
  67. }
  68. func convolveRGBASep(dst *image.RGBA, src image.Image, k *SeparableKernel) error {
  69. if len(k.X) != len(k.Y) {
  70. return fmt.Errorf("graphics: kernel not square (x %d, y %d)", len(k.X), len(k.Y))
  71. }
  72. if len(k.X)%2 != 1 {
  73. return fmt.Errorf("graphics: kernel length (%d) not odd", len(k.X))
  74. }
  75. radius := (len(k.X) - 1) / 2
  76. // buf holds the result of vertically blurring src.
  77. bounds := dst.Bounds()
  78. width, height := bounds.Dx(), bounds.Dy()
  79. buf := make([]float64, width*height*4)
  80. for y := bounds.Min.Y; y < bounds.Max.Y; y++ {
  81. for x := bounds.Min.X; x < bounds.Max.X; x++ {
  82. var r, g, b, a float64
  83. // k0 is the kernel weight for the center pixel. This may be greater
  84. // than kernel[0], near the boundary of the source image, to avoid
  85. // vignetting.
  86. k0 := k.X[radius]
  87. // Add the pixels from above.
  88. for i := 1; i <= radius; i++ {
  89. f := k.Y[radius-i]
  90. if y-i < bounds.Min.Y {
  91. k0 += f
  92. } else {
  93. or, og, ob, oa := src.At(x, y-i).RGBA()
  94. r += float64(or>>8) * f
  95. g += float64(og>>8) * f
  96. b += float64(ob>>8) * f
  97. a += float64(oa>>8) * f
  98. }
  99. }
  100. // Add the pixels from below.
  101. for i := 1; i <= radius; i++ {
  102. f := k.Y[radius+i]
  103. if y+i >= bounds.Max.Y {
  104. k0 += f
  105. } else {
  106. or, og, ob, oa := src.At(x, y+i).RGBA()
  107. r += float64(or>>8) * f
  108. g += float64(og>>8) * f
  109. b += float64(ob>>8) * f
  110. a += float64(oa>>8) * f
  111. }
  112. }
  113. // Add the central pixel.
  114. or, og, ob, oa := src.At(x, y).RGBA()
  115. r += float64(or>>8) * k0
  116. g += float64(og>>8) * k0
  117. b += float64(ob>>8) * k0
  118. a += float64(oa>>8) * k0
  119. // Write to buf.
  120. o := (y-bounds.Min.Y)*width*4 + (x-bounds.Min.X)*4
  121. buf[o+0] = r
  122. buf[o+1] = g
  123. buf[o+2] = b
  124. buf[o+3] = a
  125. }
  126. }
  127. // dst holds the result of horizontally blurring buf.
  128. for y := 0; y < height; y++ {
  129. for x := 0; x < width; x++ {
  130. var r, g, b, a float64
  131. k0, off := k.X[radius], y*width*4+x*4
  132. // Add the pixels from the left.
  133. for i := 1; i <= radius; i++ {
  134. f := k.X[radius-i]
  135. if x-i < 0 {
  136. k0 += f
  137. } else {
  138. o := off - i*4
  139. r += buf[o+0] * f
  140. g += buf[o+1] * f
  141. b += buf[o+2] * f
  142. a += buf[o+3] * f
  143. }
  144. }
  145. // Add the pixels from the right.
  146. for i := 1; i <= radius; i++ {
  147. f := k.X[radius+i]
  148. if x+i >= width {
  149. k0 += f
  150. } else {
  151. o := off + i*4
  152. r += buf[o+0] * f
  153. g += buf[o+1] * f
  154. b += buf[o+2] * f
  155. a += buf[o+3] * f
  156. }
  157. }
  158. // Add the central pixel.
  159. r += buf[off+0] * k0
  160. g += buf[off+1] * k0
  161. b += buf[off+2] * k0
  162. a += buf[off+3] * k0
  163. // Write to dst, clamping to the range [0, 255].
  164. dstOff := (y-dst.Rect.Min.Y)*dst.Stride + (x-dst.Rect.Min.X)*4
  165. dst.Pix[dstOff+0] = uint8(clamp(r+0.5, 0, 255))
  166. dst.Pix[dstOff+1] = uint8(clamp(g+0.5, 0, 255))
  167. dst.Pix[dstOff+2] = uint8(clamp(b+0.5, 0, 255))
  168. dst.Pix[dstOff+3] = uint8(clamp(a+0.5, 0, 255))
  169. }
  170. }
  171. return nil
  172. }
  173. func convolveRGBA(dst *image.RGBA, src image.Image, k Kernel) error {
  174. b := dst.Bounds()
  175. bs := src.Bounds()
  176. w := k.Weights()
  177. size, err := kernelSize(w)
  178. if err != nil {
  179. return err
  180. }
  181. radius := (size - 1) / 2
  182. for y := b.Min.Y; y < b.Max.Y; y++ {
  183. for x := b.Min.X; x < b.Max.X; x++ {
  184. if !image.Pt(x, y).In(bs) {
  185. continue
  186. }
  187. var r, g, b, a, adj float64
  188. for cy := y - radius; cy <= y+radius; cy++ {
  189. for cx := x - radius; cx <= x+radius; cx++ {
  190. factor := w[(cy-y+radius)*size+cx-x+radius]
  191. if !image.Pt(cx, cy).In(bs) {
  192. adj += factor
  193. } else {
  194. sr, sg, sb, sa := src.At(cx, cy).RGBA()
  195. r += float64(sr>>8) * factor
  196. g += float64(sg>>8) * factor
  197. b += float64(sb>>8) * factor
  198. a += float64(sa>>8) * factor
  199. }
  200. }
  201. }
  202. if adj != 0 {
  203. sr, sg, sb, sa := src.At(x, y).RGBA()
  204. r += float64(sr>>8) * adj
  205. g += float64(sg>>8) * adj
  206. b += float64(sb>>8) * adj
  207. a += float64(sa>>8) * adj
  208. }
  209. off := (y-dst.Rect.Min.Y)*dst.Stride + (x-dst.Rect.Min.X)*4
  210. dst.Pix[off+0] = uint8(clamp(r+0.5, 0, 0xff))
  211. dst.Pix[off+1] = uint8(clamp(g+0.5, 0, 0xff))
  212. dst.Pix[off+2] = uint8(clamp(b+0.5, 0, 0xff))
  213. dst.Pix[off+3] = uint8(clamp(a+0.5, 0, 0xff))
  214. }
  215. }
  216. return nil
  217. }
  218. // Convolve produces dst by applying the convolution kernel k to src.
  219. func Convolve(dst draw.Image, src image.Image, k Kernel) (err error) {
  220. if dst == nil || src == nil || k == nil {
  221. return nil
  222. }
  223. b := dst.Bounds()
  224. dstRgba, ok := dst.(*image.RGBA)
  225. if !ok {
  226. dstRgba = image.NewRGBA(b)
  227. }
  228. switch k := k.(type) {
  229. case *SeparableKernel:
  230. err = convolveRGBASep(dstRgba, src, k)
  231. default:
  232. err = convolveRGBA(dstRgba, src, k)
  233. }
  234. if err != nil {
  235. return err
  236. }
  237. if !ok {
  238. draw.Draw(dst, b, dstRgba, b.Min, draw.Src)
  239. }
  240. return nil
  241. }