b2419c92e0f01a33ec3e5d4af8deae56317039a2
[packages/hpc.git] / tests / raytrace / Illumination.hs
1 -- Copyright (c) 2000 Galois Connections, Inc.
2 -- All rights reserved. This software is distributed as
3 -- free software under the license in the file "LICENSE",
4 -- which is included in the distribution.
5
6 -- Modified to use stdout (for testing)
7
8 module Illumination
9 ( Object
10 , Light (..)
11 , light, pointlight, spotlight
12 , render
13 ) where
14
15 import Array
16 import Char(chr)
17 import Maybe
18
19 import Geometry
20 import CSG
21 import Surface
22 import Misc
23
24 type Object = CSG (SurfaceFn Color Double)
25
26 data Cxt = Cxt {ambient::Color, lights::[Light], object::Object, depth::Int}
27 deriving Show
28
29 render :: (Matrix,Matrix) -> Color -> [Light] -> Object -> Int ->
30 Radian -> Int -> Int -> String -> IO ()
31 render (m,m') amb ls obj dep fov wid ht file
32 = do { debugging
33 ; txt <- readFile "galois.sample"
34 ; let vals = read txt
35 ; let rt_vals = showBitmap' wid ht pixels
36 ; if length vals /= length rt_vals
37 then print ("BAD LENGTH",length vals,length rt_vals)
38 else do {
39 ; let cmp = sum(zipWith (\ a b -> abs (a - b) * abs (a - b)) vals rt_vals)
40 ; print $ if cmp <= (length vals * 16) then ("GOOD MATCH") else ("BAD MATCH:" ++ show cmp)
41 }}
42
43 where
44 debugging = return ()
45 {-
46 do { putStrLn (show cxt)
47 ; putStrLn (show (width, delta, aspect, left, top))
48 }
49 -}
50 obj' = transform (m',m) obj
51 ls' = [ transformLight m' l | l <- ls ]
52 pixelA = listArray ((1,1), (ht,wid))
53 [ illumination cxt (start,pixel i j)
54 | j <- take ht [0.5..]
55 , i <- take wid [0.5..] ]
56 antiA = pixelA //
57 [ (ix, superSample ix (pixelA ! ix))
58 | j <- [2 .. ht - 1], i <- [2 .. wid - 1]
59 , let ix = (j, i)
60 , contrast ix pixelA ]
61 pixels = [ [ illumination cxt (start,pixel i j) | i<- take wid [0.5..] ]
62 | j <- take ht [0.5..]
63 ]
64 cxt = Cxt {ambient=amb, lights=ls', object=obj', depth=dep}
65 start = point 0 0 (-1)
66 width = 2 * tan (fov/2)
67 delta = width / fromIntegral wid
68 aspect = fromIntegral ht / fromIntegral wid
69 left = - width / 2
70 top = - left * aspect
71 pixel i j = vector (left + i*delta) (top - j*delta) 1
72
73 superSample (y, x) col = avg $ col:
74 [ illumination cxt (start, pixel (fromIntegral x - 0.5 + xd) (fromIntegral y - 0.5 + yd))
75 | (xd, yd) <- [(-0.333, 0.0), (0.333, 0.0), (0.0, -0.333), (0.0, 0.333)]
76 ]
77
78 avg cs = divN (fromIntegral (length cs)) (uncolor (sumCC cs))
79 where divN n (r,g,b) = color (r / n) (g / n) (b / n)
80
81 contrast :: (Int, Int) -> Array (Int, Int) Color -> Bool
82 contrast (x, y) arr = any diffMax [ subCC cur (arr ! (x + xd, y + yd))
83 | xd <- [-1, 1], yd <- [-1, 1]
84 ]
85 where cur = arr ! (x, y)
86 diffMax col = (abs r) > 0.25 || (abs g) > 0.2 || (abs b) > 0.4
87 where
88 (r,g,b) = uncolor col
89
90
91 illumination :: Cxt -> Ray -> Color
92 illumination cxt (r,v)
93 | depth cxt <= 0 = black
94 | otherwise = case castRay (r,v) (object cxt) of
95 Nothing -> black
96 Just info -> illum (cxt{depth=(depth cxt)-1}) info v
97
98 illum :: Cxt -> (Point,Vector,Properties Color Double) -> Vector -> Color
99 illum cxt (pos,normV,(col,kd,ks,n)) v
100 = ambTerm `addCC` difTerm `addCC` spcTerm `addCC` recTerm
101 where
102 visibleLights = unobscured pos (object cxt) (lights cxt) normV
103 d = depth cxt
104 amb = ambient cxt
105 newV = subVV v (multSV (2 * dot normV v) normV)
106
107 ambTerm = multSC kd (multCC amb col)
108 difTerm = multSC kd (sumCC [multSC (dot normV lj) (multCC intensity col)
109 |(loc,intensity) <- visibleLights,
110 let lj = normalize ({- pos `subVV` -} loc)])
111 -- ZZ might want to avoid the phong, when you can...
112 spcTerm = multSC ks (sumCC [multSC ((dot normV hj) ** n ) (multCC intensity col)
113 |(loc,intensity) <- visibleLights,
114 -- ZZ note this is specific to the light at infinity
115 let lj = {- pos `subVV` -} normalize loc,
116 let hj = normalize (lj `subVV` normalize v)])
117 recTerm = if recCoeff `nearC` black then black else multCC recCoeff recRay
118 recCoeff = multSC ks col
119 recRay = illumination cxt (pos,newV)
120
121 showBitmapA :: Int -> Int -> Array (Int, Int) Color -> String
122 showBitmapA wid ht arr
123 = header ++ concatMap scaleColor (elems arr)
124 where
125 scaleColor col = [scalePixel r, scalePixel g, scalePixel b]
126 where (r,g,b) = uncolor col
127 header = "P6\n#Galois\n" ++ show wid ++ " " ++ show ht ++ "\n255\n"
128
129 showBitmap :: Int -> Int ->[[Color]] -> String
130 showBitmap wid ht pss
131 -- type of assert | length pss == ht && all (\ ps -> length ps == wid) pss
132 = header ++ concat [[scalePixel r,scalePixel g,scalePixel b]
133 | ps <- pss, (r,g,b) <- map uncolor ps]
134 where
135 header = "P6\n#Galois\n" ++ show wid ++ " " ++ show ht ++ "\n255\n"
136 showBitmap _ _ _ = error "incorrect length of bitmap string"
137
138 scalePixel :: Double -> Char
139 scalePixel p = chr (floor (clampf p * 255))
140
141 showBitmap' :: Int -> Int ->[[Color]] -> [Int]
142 showBitmap' wid ht pss
143 -- type of assert | length pss == ht && all (\ ps -> length ps == wid) pss
144 = concat [ concat [ [scalePixel' r,scalePixel' g,scalePixel' b]
145 | (r,g,b) <- map uncolor ps]
146 | ps <- pss ]
147 where
148 header = "P3\n#Galois\n" ++ show wid ++ " " ++ show ht ++ "\n255\n"
149 showBitmap' _ _ _ = error "incorrect length of bitmap string"
150
151 scalePixel' :: Double -> Int
152 scalePixel' p = floor (clampf p * 255)
153
154 -- Lights
155
156 data Light = Light Vector Color
157 | PointLight Point Color
158 | SpotLight Point Point Color Radian Double
159 deriving Show
160
161 light :: Coords -> Color -> Light
162 light (x,y,z) color =
163 Light (normalize (vector (-x) (-y) (-z))) color
164 pointlight (x,y,z) color =
165 PointLight (point x y z) color
166 spotlight (x,y,z) (p,q,r) col cutoff exp =
167 SpotLight (point x y z) (point p q r) col cutoff exp
168
169 transformLight m (Light v c) = Light (multMV m v) c
170 transformLight m (PointLight p c) = PointLight (multMP m p) c
171 transformLight m (SpotLight p q c r d) = SpotLight (multMP m p) (multMP m q) c r d
172
173 unobscured :: Point -> Object -> [Light] -> Vector -> [(Vector,Color)]
174 unobscured pos obj lights normV = catMaybes (map (unobscure pos obj normV) lights)
175
176 unobscure :: Point -> Object -> Vector -> Light -> Maybe (Vector,Color)
177 unobscure pos obj normV (Light vec color)
178 -- ZZ probably want to make this faster
179 | vec `dot` normV < 0 = Nothing
180 | intersects (pos `addPV` (0.0001 `multSV` vec),vec) obj = Nothing
181 | otherwise = Just (vec,color)
182 unobscure pos obj normV (PointLight pp color)
183 | vec `dot` normV < 0 = Nothing
184 | intersectWithin (pos `addPV` (0.0001 `multSV` (normalize vec)), vec) obj = Nothing
185 | otherwise = Just (vec,is)
186 where vec = pp `subPP` pos
187 is = attenuate vec color
188 unobscure org obj normV (SpotLight pos at color cutoff exp)
189 | vec `dot` normV < 0 = Nothing
190 | intersectWithin (org `addPV` (0.0001 `multSV` (normalize vec)), vec) obj = Nothing
191 | angle > cutoff = Nothing
192 | otherwise = Just (vec, is)
193 where vec = pos `subPP` org
194 vec' = pos `subPP` at
195 angle = acos (normalize vec `dot` (normalize vec'))
196
197 asp = normalize (at `subPP` pos)
198 qsp = normalize (org `subPP` pos)
199 is = attenuate vec (((asp `dot` qsp) ** exp) `multSC` color)
200
201 attenuate :: Vector -> Color -> Color
202 attenuate vec color = (100 / (99 + sq (norm vec))) `multSC` color
203
204 --
205
206 castRay ray p
207 = case intersectRayWithObject ray p of
208 (True, _, _) -> Nothing -- eye is inside
209 (False, [], _) -> Nothing -- eye is inside
210 (False, (0, b, _) : _, _) -> Nothing -- eye is inside
211 (False, (i, False, _) : _, _) -> Nothing -- eye is inside
212 (False, (t, b, (s, p0)) : _, _) ->
213 let (v, prop) = surface s p0 in
214 Just (offsetToPoint ray t, v, prop)
215
216 intersects ray p
217 = case intersectRayWithObject ray p of
218 (True, _, _) -> False
219 (False, [], _) -> False
220 (False, (0, b, _) : _, _) -> False
221 (False, (i, False, _) : _, _) -> False
222 (False, (i, b, _) : _, _) -> True
223
224 intersectWithin :: Ray -> Object -> Bool
225 intersectWithin ray p
226 = case intersectRayWithObject ray p of
227 (True, _, _) -> False -- eye is inside
228 (False, [], _) -> False -- eye is inside
229 (False, (0, b, _) : _, _) -> False -- eye is inside
230 (False, (i, False, _) : _, _) -> False -- eye is inside
231 (False, (t, b, _) : _, _) -> t < 1.0