Abstract:
The present invention provides an artifact-free, real-time global illumination rendering system. 
     Artifact-free, real-time global illumination rendering can be performed by estimating appropriate specular roughness using sampling results and performing shading using the estimated specular roughness.

Description:
BACKGROUND 
       [0001]    1. Technical Field 
         [0002]    The present invention relates to a real-time global illumination rendering system in which artifacts are suppressed by estimating appropriate specular roughness using sampling results and performing shading using the estimated specular roughness. 
         [0003]    2. Related Art 
         [0004]    In computer graphics, in order to determine brightness of an object close to the actual one, rendering by global illumination (GI) is performed. The global illumination is an illumination calculation taking into account not only a direct light component from a light source, but also an indirect light component. To represent in which direction and at what intensity the indirect light is reflected, a bidirectional reflectance distribution function (BRDF) is used. 
         [0005]    In global illumination, as a method for determining an indirect light component, there is known final gathering (HACHISUKA, T. 2005. High-quality global illumination rendering using rasterization. In GPU Gems 2. Addison-Wesley Professional, ch. 38, 615-634). Final gathering is a method for determining an indirect light component by tracing a plurality of rays. 
         [0006]    In addition, as a technique for determining an indirect light component, there is known an instant radiosity method using a virtual point light source (KELLER, A. 1997. Instant radiosity. In ACM SIGGRAPH &#39;97, 49-56). As a technique for tracing rays in the same direction at a time, there are ray-bundles. 
         [0007]    To represent in which direction and at what intensity light is reflected, a bidirectional reflectance distribution function (BRDF) is proposed. Examples of already-known BRDF models include a Lambert model representing diffuse reflection, a Phong model representing specular reflection, a Blinn model, and a Cook-Torrance model. 
       SUMMARY OF THE INVENTION 
       [0008]    To attempt to implement rendering by global illumination (GI) based on the conventional methods, there is a problem of the occurrence of artifacts. In addition, ray-bundles have a problem of poor sampling efficiency for specular reflection. In the case of off-line rendering, by increasing the number of samples, computer graphics close to the actual one can be obtained even by ray-bundles. However, since the number of calculations increases, ray-bundles are not suitable for real-time rendering. 
         [0009]    When, for a specular reflection model, shading is performed using ray-bundles or instant radiosity, there is a problem of the occurrence of sharp artifacts. 
         [0010]    An object of the present invention is therefore to provide an artifact-free, real-time global illumination rendering system. 
         [0011]    The present invention is basically based on the finding that artifact-free, real-time global illumination rendering can be performed by estimating appropriate specular roughness using sampling results and performing shading using the estimated specular roughness (α i ). 
         [0012]    A first aspect of the present invention is directed to a shading method by a computer, in which appropriate specular roughness that does not cause artifacts is estimated using sampling results and shading is performed using the estimated specular roughness. This method includes a sampling step (S 101 ) and an estimated specular roughness obtaining step (S 102 ). After obtaining estimated specular roughness, shading is performed based on a normal method. 
         [0013]    The sampling step (S 101 ) is a step of sampling, for a plurality of samples, a path of light reflected at a given object from a virtual light source and reaching a given eye point. 
         [0014]    The estimated specular roughness obtaining step (S 102 ) is a step for obtaining estimated specular roughness of a surface of the object by obtaining a distribution function for the surface of the object using the paths obtained in the sampling step and performing arithmetic processing such that a maximum value of the distribution function is less than or equal to a predetermined value. By using the estimated specular roughness obtained in this step, artifact-free, real-time global illumination rendering can be performed. 
         [0015]    A second aspect of the present invention is directed to a computer graphics system that estimates, using sampling results, appropriate specular roughness that does not cause artifacts, and performs shading. That is, the system is a computer graphics system that implements the sampling method of the first aspect. The system includes a sampling unit  101  and an estimated specular roughness obtaining unit  102 . 
         [0016]    The sampling unit  101  is a component for sampling, for a plurality of samples, a path of light reflected at a given object from a virtual light source and reaching a given eye point. 
         [0017]    The estimated specular roughness obtaining unit  102  is a component for obtaining estimated specular roughness of a surface of the object by performing arithmetic processing using the paths obtained by the sampling unit such that a maximum value of a distribution function is less than or equal to a predetermined value. 
         [0018]    A shading method of the present invention is based on the above-described method and thus can determine estimated specular roughness (α i ) of a surface of an object with suppressed artifacts based on a small number of sampling data and can perform shading. Thus, with the use of the shading method of the present invention, real-time global illumination rendering can be implemented in particular for scenes including specular reflection. 
     
    
     
       BRIEF DESCRIPTION OF DRAWINGS 
         [0019]      FIG. 1  is a flowchart for describing an example of a shading method related to a first embodiment of the present invention; 
           [0020]      FIG. 2  is a block diagram of a computer graphics system according to the first embodiment; 
           [0021]      FIG. 3  is a diagram for describing a light path of indirect illumination; 
           [0022]      FIG. 4  is a block diagram showing an exemplary hardware configuration of a system of the present invention; 
           [0023]      FIG. 5  is a flowchart showing an example of steps for describing a shading method related to a second embodiment of the present invention; 
           [0024]      FIG. 6  is a block diagram for describing a computer graphics rendering system using a computer which is related to the second embodiment of the present invention; and 
           [0025]      FIGS. 7A to 7D  are diagrams showing a point light source and specular image computer graphics in a Cornell box. 
       
    
    
     DETAILED DESCRIPTION 
       [0026]      FIG. 1  is a flowchart for describing an example of a shading method related to a first embodiment of the present invention. As shown in  FIG. 1 , the method includes a sampling step (S 101 ) and an estimated specular roughness obtaining step (S 102 ). Shading using a specular roughness parameter (α) is already known. Therefore, in the present invention, after obtaining estimated specular roughness (α i ), shading is performed based on a normal method, using the estimated specular roughness (α i ) as a specular roughness parameter. 
         [0027]      FIG. 2  is a block diagram of a computer graphics system according to the first embodiment. The computer graphics system shown in  FIG. 2  is a system that estimates, using sampling results, appropriate specular roughness that does not cause artifacts, and performs shading. The system shown in  FIG. 2  includes a sampling unit  101  and an estimated specular roughness obtaining unit  102 . Though not particularly shown, the system preferably includes a display unit and a shading unit that performs shading using obtained estimated specular roughness (α i ). By this, the system can obtain computer graphics using the obtained estimated specular roughness (α i ), and display the computer graphics on the display unit. 
         [0028]      FIG. 3  is a diagram for describing a light path of indirect illumination. Light from a light source  11 , for example, is reflected at a given point  15  on a wall  13  and is then reflected at a given point  19  on a surface of an object  17 . The point  15  serves as a virtual light source. Then, the light reflected at the given point  19  reaches an eye point  20 . This light path can be determined by, for example, a computer reading the coordinates of the eye point  20 , the surface of the object  17 , the wall  13 , and the light source  11  from a storage unit and performing a computation for determining a light path. This system is configured, for example, to be able to design coefficients for reflection in association with coordinate values. Hence, the computer stores, for example, in the storage unit, coefficients for reflection associated with the given point  19  on the object  17  and the given point  15  on the wall  13 . The coefficients, for example, may be stored in a database or may be temporarily stored in a buffer. Then, the computer reads the coefficients for reflection for the given point  19  on the object  17  and the given point  15  on the wall  13 . The computer may be configured to be able to determine reflection paths using the read coefficients. In this manner, the computer can sample a path of light reflected at the given object  17  from the virtual light source  15  and reaching the given eye point 
         [0029]    Note that upon actual calculation, a light path may be determined where a plurality of rays are emitted from the eye point  20  and light reflected at the given point  19  is reflected at the given point  15  on the wall  13  and then reaches the light source  11 . In addition, a virtual screen may be determined between the eye point  20  and the object  17 . Then, samples on the screen may be connected to the eye point  20 , and a plurality of lines of sight may be determined thereby. Furthermore, of the plurality of lines of sight, a line of sight reaching the object  17  may be determined. A method for determining the light path is already known in, for example, a ray tracing method. In addition, in actual calculation, light paths may be determined where a plurality of rays are emitted from the eye point  20  and light reflected at a plurality of locations is reflected at other points and then reaches the light source  11 . 
         [0030]    The sampling unit  101  is a component for sampling, for a plurality of samples (sample points), a path of light reflected at the given object  17  from the light source  11  and reaching the given eye point  20 . The component is implemented by a program, a control apparatus, an arithmetic apparatus, a storage apparatus, etc. A sample may be, for example, a point on a virtual screen in a path of line of sight. 
         [0031]    The estimated specular roughness obtaining unit  102  is a component for obtaining estimated specular roughness (α i ) of a surface of the object  17 , using the paths obtained by the sampling unit  101 . The estimated specular roughness obtaining unit  102  determines an intersection point (collision point) of a line of sight and the object using a path obtained by the sampling unit  101 , and determines a distribution function for the intersection point. The distribution function is stored in the storage unit so as to be associated with coordinates. Therefore, the distribution function can be read from the storage unit using the coordinates of the intersection point. The distribution function is represented by a predetermined function and is represented using an object surface roughness parameter α. The distribution function can be said to be a function indicating the roughness of the surface of the object. The object surface roughness parameter α is known and is, for example, stored in the storage unit so as to be associated with coordinate values or the object. Then, the estimated specular roughness obtaining unit  102  determines density of the paths and represents the distribution function using the density of the paths and performs arithmetic processing such that the maximum value of the distribution function is less than or equal to a predetermined value. The density of paths may be obtained analytically. In addition, the density of paths may be, for example, the density of paths present within a predetermined range from a given sample point, or maybe the density of paths present within a predetermined range from a given point on a virtual screen. In addition, the density of paths can also be determined by determining the number of paths passing through a screen with a certain area and using the area of the screen and the number of paths. To represent the distribution function using the density of paths, a relational expression between the distribution function and the density of paths is stored, and a computation for representing the distribution function is allowed to be performed based on the relational expression after determining the density of paths. This computation can be implemented by a program or can also be implemented by hardware. In this manner, the estimated specular roughness obtaining unit  102  can obtain the estimated specular roughness (α i ) of the surface of the object. The predetermined value is read from the storage unit before performing arithmetic processing. The density of paths is thus determined, and estimated specular roughness (α i ) at which the maximum value of the distribution function is less than or equal to the predetermined value is determined, and then shading is performed using the estimated specular roughness (α i ). Hence, the system of the present invention can suppress artifacts. The estimated specular roughness obtaining unit  102  is actually implemented by a program, a control apparatus, an arithmetic apparatus, a storage apparatus, etc. 
         [0032]      FIG. 4  is a block diagram showing an exemplary hardware configuration of a system of the present invention. The system includes an input apparatus  21 , a control apparatus  22 , an arithmetic apparatus  23 , a storage apparatus  24 , and an image processing block  25 . The components are connected to each other by a bus  26  so that they can give and receive information to/from each other. In addition, the system is connected to an interface (I/F)  27  via the bus  26 . Hence, for example, the system can be connected to an information recording medium  28  that stores a program, through the I/F  27 . The program is to cause a computer to function as an apparatus having predetermined means and predetermined functions. In addition, the program is to cause the computer to perform predetermined steps. Note that a memory connected through the I/F  27  may function as all or part of the storage apparatus  24 . 
         [0033]    The image processing block  25  includes a graphics processing unit (GPU)  29  and a video RAM (VRAM)  30 . The GPU  29  and the VRAM  30  are connected to each other so that they can give and receive information to/from each other. Reference numeral  31  in the drawing indicates a display screen (monitor). The GPU  29  and the VRAM  30  may implement a shading method of the present invention. 
         [0034]    In addition, for example, the working space of the control apparatus  22 , the arithmetic apparatus  23 , and the storage apparatus  24  and the program stored in the information recording medium  28  may function as the sampling unit  101  and the estimated specular roughness obtaining unit  102 . In this case, the storage areas of the storage apparatus  24  and the information recording medium  28  function as areas that store information necessary for computations and working areas. When operation information is inputted from the input apparatus  21 , the operation information is passed to the control apparatus  22  via the bus  26 . The control apparatus  22  then reads the program stored in the information recording medium  28  and performs a predetermined process. The control apparatus  22  reads various information stored in the storage apparatus  24  and the information recording medium  28  in response to an instruction from the program, and the arithmetic apparatus  23  performs a predetermined computation. The control apparatus  22  appropriately stores the results of the computation in the storage apparatus  24 , performs a predetermined process using the results of the computation, and appropriately stores the resulting information in the storage apparatus  24  such as a frame buffer, and then, outputs appropriate information from the monitor  31 . 
         [0035]    The above-described system may function as a geometry unit, a rendering unit, a frame buffer, and a display unit in order to implement a shading method of the present invention to obtain computer graphics. The geometry unit receives polygon data and performs calculations including a brightness calculation. The rendering unit performs a rendering process using vertex information and information including brightness data. Raster data obtained through the rendering unit is temporarily stored in the frame buffer and is displayed on the display unit. 
         [0036]    Next, a first sampling method of the present invention will be described. 
         [0037]    A sampling step (S 101 ) is a step of sampling, for a plurality of samples, a path of light reflected at the given object  17  from the virtual light source  15  and reaching the given eye point  20 . For the sampling step, for example, a normal sampling step in a ray tracing method can be adopted. 
         [0038]    An example of the sampling step will be described below. Note that the following example is merely exemplification and thus it is apparent that any other method than the following method may be used as long as a path that reaches the eye point  20  from the light source  11  through the virtual light source  15  and the object  17  can be determined. 
         [0039]    Step  1 : Lines of sight going in a direction of the object  17  from the given eye point  20  through a given pixel on a screen are determined. Thereafter, the object  17  intersecting the determined lines of sight is checked. In this step, for example, the coordinates of the eye point  20  are read from the storage unit. Then, from the coordinates of the eye point  20 , lines of sight are radiated in a plurality of directions. At this time, lines of sight passing through a given pixel on the screen may be determined. An object(s) intersecting the respective lines of sight is(are) checked. At this time, the coordinates of one or a plurality of objects (or the vertex coordinates of polygons composing an object(s)) are read from the storage unit. Then, using the read coordinates, a computation for determining whether the lines of sight intersect the object(s) is performed. In this manner, it can be determined whether lines of sight passing through a given pixel on the screen intersect a given object. After determining whether the given object intersects the lines of sight, information on the intersecting object is stored in the storage unit for each line of sight. 
         [0040]    Step  2 : When there is an object intersecting a given line of sight, a computation for determining an intersection point of the line of sight and the object is performed. When there are a plurality of objects intersecting a given line of sight, for example, intersection points of the given line of sight and all of the objects maybe determined. At this time, the vertex coordinates of polygons composing an object may be read from the storage unit, and a computation for determining an intersection point of a given polygon and a given line of sight may be performed. 
         [0041]    Step  3 : Distances to intersection points of a given line of sight and objects intersecting the given line of sight may be determined, and an object closest to the eye point may be extracted. Then, information indicating that, of the objects intersecting the given line of sight, a given object is present closest to the eye point side may also be stored together. The fact that a given object is present closest to the eye point side for a given line of sight can be determined by, for example, performing a computation for comparing distances between the eye point and the coordinates of the intersection points of the line of sight and objects. 
         [0042]    In this manner, a group of lines of sight intersecting a given object can be determined. Then, information about those of the plurality of lines of sight that are determined to intersect the given object may be stored in the storage unit. Furthermore, information about a given object that intersects one of the plurality of lines of sight and that is, moreover, determined to be present at the forefront may be stored in the storage unit. By reading these pieces of information stored in the storage unit, the lines of sight intersecting the given object can be read. 
         [0043]    It is determined whether an object or a wall present in the above-described path is an object causing reflection or refraction. Since the storage unit stores information about objects (or polygons included in objects) and information about whether the objects are objects causing reflection or objects causing refraction, such a determination is made by reading the information from the storage unit. When the object or wall present in the above-described path is an object causing reflection or refraction, data for determining a reflection direction or a refraction direction is read from the storage unit, and a computation is performed using the read data together with information about the above-described direction of the line of sight, by which a reflection direction or a refraction direction is determined. Then, using the determined reflection direction or refraction direction as a new line of sight, the above-described steps  2  and  3  are performed and a computation for extracting an object that appears reflected or refracted is performed. 
         [0044]    By doing the above, a plurality of paths on a given screen can be determined. In addition, a path for one of a plurality of samples on a given screen that is related to a line of sight intersecting a given object can also be determined. 
         [0045]    The estimated specular roughness obtaining step (S 102 ) is a step for obtaining estimated specular roughness (α i ) of the surface of the object by obtaining a distribution function for the surface of the object using the plurality of paths obtained in the sampling step, and performing arithmetic processing such that the maximum value of the distribution function is less than or equal to a predetermined value. 
         [0046]    For example, in this step, the computer reads, from the storage unit, paths for a group of lines of sight intersecting a given object among lines of sight from a given eye point, which are determined in the sampling step. Then, a distribution function for each intersection point associated with intersection point coordinates is read from the storage unit. The distribution function is a known function indicating the reflectance distribution of incident light. 
         [0047]    When an intersection point (collision point) distribution function is determined, estimated specular roughness (α i ) of the surface of the object is obtained. Specifically, from the plurality of paths determined in the sampling step, the density of the paths is determined and a distribution function is represented. Then, arithmetic processing is performed such that the maximum value of the distribution function is less than or equal to a predetermined value represented including the density of the paths. To implement this step, for example, a predetermined value is stored in the storage unit in advance. When this computation is performed, the predetermined value is read from the storage unit. Then, a computation for determining the maximum value of the distribution function is performed. Thereafter, a computation is performed such that the maximum value of the distribution function is less than or equal to (or is less than) the predetermined value. By doing so, estimated specular roughness (α i ) at which the maximum value of the distribution function is less than or equal to the predetermined value can be determined. The system of the present invention performs shading using the estimated specular roughness (α i ). Hence, a situation in which artifacts occur can be suppressed. 
         [0048]    Next, a second embodiment of the present invention will be described. The second embodiment is an embodiment for definitely achieving determination of appropriate estimated specular reflection roughness that does not cause artifacts, by a shading method of the present invention.  FIG. 5  shows an example of steps for describing a shading method related to the second embodiment of the present invention. As shown in  FIG. 5 , the shading method includes a sampling step (S 101 ), an estimated density calculation step (S 112 ), a distribution function obtaining step (S 113 ), and an estimated specular roughness calculation step (S 114 ). 
         [0049]      FIG. 6  is a block diagram for describing a computer graphics rendering system using a computer which is related to the second embodiment of the present invention. As shown in  FIG. 6 , the system includes a sampling unit  101 , an estimated density calculating unit  112 , a distribution function obtaining unit  113 , and an estimated specular roughness calculating unit  114 . The sampling unit  101 , the estimated density calculating unit  112 , the distribution function obtaining unit  113 , and the estimated specular roughness calculating unit  114  are virtual components of a computer for implementing a sampling step (S 101 ), an estimated density calculation step (S 112 ), a distribution function obtaining step (S 113 ), and an estimated specular roughness calculation step (S 114 ), respectively. The computer&#39;s input/output unit, storage unit, arithmetic unit, control unit, and bus connecting the units implement the above-described components. 
         [0050]    The sampling unit  101  is a means for sampling light paths for a plurality of samples reflected at an object from a virtual light source and reaching a given eye point. 
         [0051]    The estimated density calculating unit  112  is a means for determining, using the paths obtained by the sampling unit  101 , estimated density (probability density) which is the density of the paths. 
         [0052]    The distribution function obtaining unit  113  is a means for obtaining a distribution function for an intersection point of the object and one of the paths determined by the sampling unit  101 . 
         [0053]    The estimated specular roughness calculating unit  114  is a means for determining estimated specular roughness of the surface of the object by performing arithmetic processing on the distribution function obtained by the distribution function obtaining unit  113 , using the density of the paths determined by the estimated density calculating unit  112  such that the maximum value is less than or equal to a predetermined value. 
         [0054]    The sampling step (S 101 ) is a step for sampling, by the computer, light paths for samples reflected at a given object from a virtual light source and reaching a given eye point. The virtual light source is a virtual light source in three-dimensional computer graphics. This step is the same as the sampling step of the first embodiment. 
         [0055]    The estimated density calculation step (S 112 ) is a step for determining, using the paths obtained in the sampling step, estimated density (probability density) which is the density of the paths. The estimated density at position y can be represented as d hat(y). In the sampling step, normally, a plurality of paths are determined. In this example, the density of the paths is determined analytically. Specifically, when a plurality of paths are determined in advance, the density of the paths may be able to be determined according to the number of the determined paths, or the density of the paths may be able to be determined from the positions of the determined plurality of paths. In the former case, after determining the number of paths, predetermined arithmetic processing is performed using the determined number of paths, by which the density of the paths can be determined. In the latter case, when a plurality of paths together with their positions are determined in the sampling step, predetermined arithmetic processing is performed, by which the density of the paths can be determined. 
         [0056]    The distribution function obtaining step (S 113 ) is a step for obtaining a distribution function for an intersection point of the object and one of the paths determined in the sampling step (S 101 ). In the sampling step (S 101 ), an intersection point of a line of sight (i.e., a path) and the object is determined. Then, in this step (S 113 ) ; for example, the computer reads the coordinates of the intersection point from the storage unit. Then, the computer reads, using the coordinates of the intersection point, a distribution function stored so as to be associated with the read coordinates of the intersection point from the storage unit. In this manner, the distribution function for the intersection point of the path and the object can be obtained. Note that the distribution function is a predetermined function using the roughness parameter α of the object. 
         [0057]    The estimated specular roughness calculation step (S 114 ) is a step for determining estimated specular roughness of the surface of the object including a sample, by performing arithmetic processing such that the maximum value of the distribution function obtained in the distribution function obtaining step (S 113 ) is less than or equal to a predetermined value. The estimated specular roughness (α i ) is a value that estimates specular roughness at which artifacts are less likely to occur and that is used when determining specular reflection at a surface of a given object in computer graphics. In the present invention, by using, as surface roughness, the estimated specular roughness instead of already known surface roughness (α), computer graphics with suppressed artifacts can be achieved. The distribution function can be represented using, for example, the density of the paths d hat (y). Then, by performing arithmetic processing such that the maximum value of the distribution function is less than or equal to the predetermined value, estimated specular roughness (α i ) can be determined. For example, the maximum value of the distribution function (i.e., the maximum value of a function including the density of the paths) is determined, and a coefficient of the function including the density of the paths is determined such that the maximum value is less than or equal to a predetermined value. Thereafter, surface roughness is represented using the coefficient and the density of the paths. Then, a computation is performed to determine the smaller one of the surface roughness using the coefficient and the density of the paths and the already known surface roughness (α). Then, in this embodiment, the larger one of the parameter of the determined smaller one and  0  is determined to be estimated specular roughness. This computation can be performed by a circuit or a program that reads, by the computer, a distribution function and a density of the paths d hat (y) from the storage unit and performs a predetermined computation. 
         [0058]    The computer can determine estimated specular roughness in the above-described manner. Then, the computer performs a shading process in which a scene including indirect light is rendered using the estimated specular roughness. The system of the present invention is particularly effective for rendering a surface of a specular object. After performing the shading process, rasterization is performed, by which computer graphics can be obtained. 
       EXAMPLE 1 
       [0059]    The method of the following example compresses artifacts by smoothing. The method calculates the roughness parameter of gloss article, which has brilliance, by clamping based on assumption of kernel density. Generally, density estimation is represented by: 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       d 
                       ^ 
                     
                      
                     
                       ( 
                       y 
                       ) 
                     
                   
                   = 
                   
                     
                       ∑ 
                       i 
                       N 
                     
                      
                     
                       
                         
                           K 
                            
                           
                             ( 
                             
                               
                                  
                                 
                                    
                                   
                                     
                                       y 
                                       i 
                                     
                                     - 
                                     y 
                                   
                                    
                                 
                                  
                               
                               , 
                               r 
                             
                             ) 
                           
                         
                         
                           V 
                            
                           
                             ( 
                             r 
                             ) 
                           
                         
                       
                        
                       
                         Φ 
                          
                         
                           ( 
                           
                             y 
                             i 
                           
                           ) 
                         
                       
                     
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     1 
                   
                   ] 
                 
               
             
           
         
       
     
         [0060]    where ̂ d(y) (d hat(y)) is the estimated density at the position y, K the kernel function, y i  the position of the ith particle, Φ the weight of the particle, r the kernel radius, and V(r) is a normalization factor. 
         [0061]    The particles correspond to calculated paths. 
         [0062]    If the radius is too small, the estimated density distribution can be spiky and extremely different from the actual one. We may use a different kernel radius for each particle. 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       d 
                       ^ 
                     
                      
                     
                       ( 
                       y 
                       ) 
                     
                   
                   = 
                   
                     
                       ∑ 
                       i 
                       N 
                     
                      
                     
                       
                         
                           K 
                            
                           
                             ( 
                             
                               
                                  
                                 
                                    
                                   
                                     
                                       y 
                                       i 
                                     
                                     - 
                                     y 
                                   
                                    
                                 
                                  
                               
                               , 
                               
                                 r 
                                 i 
                               
                             
                             ) 
                           
                         
                         
                           V 
                            
                           
                             ( 
                             
                               r 
                               i 
                             
                             ) 
                           
                         
                       
                        
                       
                         Φ 
                          
                         
                           ( 
                           
                             y 
                             i 
                           
                           ) 
                         
                       
                     
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     2 
                   
                   ] 
                 
               
             
           
         
       
     
         [0063]    If the density ̂ d(y i ) is given, we can obtain the appropriate r i  using the following equation: 
         [0000]    
       
         
           
             
               
                 
                   
                     r 
                     i 
                   
                   = 
                   
                     
                       V 
                       
                         - 
                         1 
                       
                     
                      
                     
                       ( 
                       
                         n 
                         
                           
                             d 
                             ^ 
                           
                            
                           
                             ( 
                             
                               y 
                               i 
                             
                             ) 
                           
                         
                       
                       ) 
                     
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     3 
                   
                   ] 
                 
               
             
           
         
       
     
         [0064]    where n is the number of nearest neighbor particles. We treat the distribution function D as the kernel function as: 
         [0000]    
       
         
           
             
               
                 
                   
                     D 
                      
                     
                       ( 
                       θ 
                       ) 
                     
                   
                   = 
                   
                     
                       K 
                        
                       
                         ( 
                         
                           θ 
                           , 
                           r 
                         
                         ) 
                       
                     
                     
                       V 
                        
                       
                         ( 
                         r 
                         ) 
                       
                     
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     4 
                   
                   ] 
                 
               
             
           
         
       
     
         [0065]    That is, we consider the roughness parameter of a material approximately corresponds to r. The maximum value of D is limited as follows: 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       max 
                       θ 
                     
                      
                     
                       
                         ( 
                         
                           D 
                            
                           
                             ( 
                             θ 
                             ) 
                           
                         
                         ) 
                       
                        
                       
                         n 
                         
                           
                             d 
                             ^ 
                           
                            
                           
                             ( 
                             
                               y 
                               i 
                             
                             ) 
                           
                         
                       
                     
                   
                   ≤ 
                   1 
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     5 
                   
                   ] 
                 
               
             
           
         
       
     
         [0066]    To satisfy the above Equation, we provide a new concept regarding the roughness of the material. 
         [0067]    For example, we show the case of D being phong distribution as follows: 
         [0000]    
       
         
           
             
               
                 
                   
                     D 
                      
                     
                       ( 
                       θ 
                       ) 
                     
                   
                   = 
                   
                     
                       
                         α 
                         + 
                         1 
                       
                       
                         2 
                          
                         π 
                       
                     
                      
                     
                       
                         cos 
                         α 
                       
                        
                       
                         ( 
                         θ 
                         ) 
                       
                     
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     6 
                   
                   ] 
                 
               
             
           
         
       
     
         [0068]    where α (alpha) is the phong exponent to control roughness of a surface. In this case, alpha meets following: 
         [0000]    
       
         
           
             
               
                 
                   0 
                   ≤ 
                   α 
                   ≤ 
                   
                     
                       
                         2 
                          
                         π 
                          
                         
                             
                         
                          
                         
                           
                             d 
                             ^ 
                           
                            
                           
                             ( 
                             
                               y 
                               i 
                             
                             ) 
                           
                         
                       
                       n 
                     
                     - 
                     1 
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     7 
                   
                   ] 
                 
               
             
           
         
       
     
         [0069]    Then, we can clamp roughness parameter α of the ith particle as the following estimated roughness parameter α i . 
         [0000]    
       
         
           
             
               
                 
                   
                     α 
                     i 
                   
                   = 
                   
                     max 
                     ( 
                     
                       0 
                       , 
                       
                         min 
                         ( 
                         
                           α 
                           , 
                           
                             
                               
                                 2 
                                  
                                 π 
                                  
                                 
                                     
                                 
                                  
                                 
                                   
                                     d 
                                     ^ 
                                   
                                    
                                   
                                     ( 
                                     
                                       y 
                                       i 
                                     
                                     ) 
                                   
                                 
                               
                               n 
                             
                             - 
                             1 
                           
                         
                         ) 
                       
                     
                     ) 
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     8 
                   
                   ] 
                 
               
             
           
         
       
     
         [0070]    To clamp the roughness using the above distribution function D, we use the following assumption: 
         [0000]    
       
         
           
             
               
                 
                   
                     1 
                     
                       
                         d 
                         ^ 
                       
                        
                       
                         ( 
                         
                           y 
                           i 
                         
                         ) 
                       
                     
                   
                   = 
                   
                     
                       1 
                       
                         d 
                          
                         
                           ( 
                           
                             x 
                             , 
                             
                               ω 
                               i 
                             
                           
                           ) 
                         
                       
                     
                     + 
                     
                       1 
                       
                         d 
                          
                         
                           ( 
                           
                             x 
                             , 
                             
                               ω 
                               i 
                               ′ 
                             
                           
                           ) 
                         
                       
                     
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     9 
                   
                   ] 
                 
               
             
           
         
       
     
         [0071]    where d(x,w i ) is the density of samples per unit solid angle at the position x and direction w i . The w i  and w′ i  are the incoming and outgoing direction of the ith sample, respectively. We consider the situation when the path (x 0 , x 1 , x 2 , x 3 ) is sampled (see  FIG. 4 ). For the sake of simplicity, we write ̂ d(yj) at x j  as ̂ d j , and d(x j ,w j ,k) as d j,k . The wi,j means the direction from the position x j  to the position x k . 
         [0072]    In calculation of FG, final gathering, considering that D 3,4  is infinity, ̂ d is given by: 
         [0000]        {circumflex over (d)}   3   =d   3,2   =Np   g (ω 3,2 )   [Number 10]
 
         [0073]    N indicates the number of samples and p g  indicates PDF (probability density function) 
         [0074]    In the case of one-bounce indirect illumination from a point light represented by the path (x 0 , x j , x 3 ) where j={1; 2}, ̂ d is obtained as follows: 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       
                         d 
                         
                           0 
                           , 
                           j 
                         
                       
                       = 
                       
                         
                           d 
                           
                             j 
                             , 
                             0 
                           
                         
                         = 
                         
                           N 
                            
                           
                             ( 
                             
                               
                                 
                                   
                                     p 
                                     g 
                                   
                                    
                                   
                                     ( 
                                     
                                       ω 
                                       
                                         3 
                                         , 
                                         j 
                                       
                                     
                                     ) 
                                   
                                 
                                  
                                 
                                   
                                     
                                        
                                       
                                          
                                         
                                           
                                             x 
                                             j 
                                           
                                           - 
                                           
                                             x 
                                             0 
                                           
                                         
                                          
                                       
                                        
                                     
                                     2 
                                   
                                   
                                     
                                        
                                       
                                          
                                         
                                           
                                             x 
                                             j 
                                           
                                           - 
                                           
                                             x 
                                             3 
                                           
                                         
                                          
                                       
                                        
                                     
                                     2 
                                   
                                 
                               
                               + 
                               
                                 
                                   p 
                                   e 
                                 
                                  
                                 
                                   ( 
                                   
                                     ω 
                                     
                                       0 
                                       , 
                                       j 
                                     
                                   
                                   ) 
                                 
                               
                             
                             ) 
                           
                         
                       
                     
                     ; 
                   
                    
                   
                     
 
                   
                    
                   
                     
                       d 
                       
                         3 
                         , 
                         j 
                       
                     
                     = 
                     
                       
                         d 
                         
                           j 
                           , 
                           3 
                         
                       
                       = 
                       
                         N 
                          
                         
                           ( 
                           
                             
                               
                                 p 
                                 g 
                               
                                
                               
                                 ( 
                                 
                                   ω 
                                   
                                     3 
                                     , 
                                     j 
                                   
                                 
                                 ) 
                               
                             
                             + 
                             
                               
                                 
                                   p 
                                   e 
                                 
                                  
                                 
                                   ( 
                                   
                                     ω 
                                     
                                       0 
                                       , 
                                       j 
                                     
                                   
                                   ) 
                                 
                               
                                
                               
                                 
                                   
                                      
                                     
                                        
                                       
                                         
                                           x 
                                           j 
                                         
                                         - 
                                         
                                           x 
                                           3 
                                         
                                       
                                        
                                     
                                      
                                   
                                   2 
                                 
                                 
                                   
                                      
                                     
                                        
                                       
                                         
                                           x 
                                           j 
                                         
                                         - 
                                         
                                           x 
                                           0 
                                         
                                       
                                        
                                     
                                      
                                   
                                   2 
                                 
                               
                             
                           
                           ) 
                         
                       
                     
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     11 
                   
                   ] 
                 
               
             
             
               
                 
                   
                     
                       
                         d 
                         ^ 
                       
                       j 
                     
                     = 
                     
                       1 
                       
                         
                           1 
                           
                             d 
                             
                               j 
                               , 
                               0 
                             
                           
                         
                         + 
                         
                           1 
                           
                             d 
                             
                               j 
                               , 
                               3 
                             
                           
                         
                       
                     
                   
                    
                   
                     
 
                   
                    
                   
                     
                       
                         d 
                         ^ 
                       
                       3 
                     
                     = 
                     
                       
                         d 
                         
                           3 
                           , 
                           j 
                         
                       
                       . 
                     
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     12 
                   
                   ] 
                 
               
             
           
         
       
     
         [0075]    The appropriate n can be determined from the weight of reflection. 
         [0000]    
       
         
           
             
               
                 
                   n 
                   = 
                   
                     
                       n 
                       ′ 
                     
                      
                     
                       M 
                       
                         D 
                          
                         
                           ( 
                           θ 
                           ) 
                         
                       
                     
                   
                 
               
               
                 
                   [ 
                   
                     Number 
                      
                     
                         
                     
                      
                     10 
                   
                   ] 
                 
               
             
           
         
       
     
         [0076]    where n′ is the user-specified parameter. 
       EXAMPLE  2   
       [0077]    The above-described algorithm is evaluated. The results of the evaluation are shown in  FIGS. 7A to 7D .  FIGS. 7A to 7D  are diagrams showing a point light source and specular image computer graphics in a Cornell box (α=300 and N=84).  FIG. 7A  shows an example in which suppression is not performed (comparative example),  FIG. 7B  shows an example of VSLs (maximum 10 sample points per VSL) (comparative example),  FIG. 7C  shows exemplary rendering of the present invention (parameter n′ related to the maximum number of virtual neighbors n=8), and  FIG. 7D  shows the ground truth (using, as the frame buffer, 640×480 pixels, VPL shadow map: 256×256 pixels, and ray-bundles: 256×256 pixels and using as, the GPU, an AMD RADEON HD6850). 
         [0078]    Since the same random numbers are used for all pixels in this implementation example, too, some artifacts occur due to the variance of Monte Carlo integration. 
         [0079]    The present invention relates to computer graphics and thus can be used in computer and amusement fields. 
         [0080]    Note that the present invention also includes matters included in the accompanying paper and the documents cited in the paper, and an invention obvious to those skilled in the art therefrom.