A method of adaptive direct volume rendering, comprising: fragmenting a sampled 3-D dataset of a scalar field into a plurality of sub-volumes of different sizes, each sub-volume associated with a set of data value parameters characterizing the data value distribution of the scalar field within the sub-volume; defining an opacity transfer function that specifies a relationship between data values of the scalar field and opacity level; selectively casting a plurality of rays from a 2-D image plane towards the sampled dataset (106; 302), each ray having an initial ray energy and a cross-section; for each ray cast from a selected location on the 2-D image plane, adaptively selecting a subset of the plurality of sub-volumes for interaction with the ray based on the data value distribution of the scalar field within the sub-volumes that are located along the path of the ray; estimating the portion of the ray energy reflected by each sub-volume of the subset using an illumination model, the data values in the sub-volume and the opacity specified by the opacity transfer function for the data values in the sub-volume; and summing the reflected portions of the ray energy as a pixel value at the selected location on the 2-D image plane; and estimating pixel values at other locations on the 2-D image plane using the pixel values at the selected locations, wherein the step of selecting a subset of the plurality of sub-volumes for interacting with the ray includes: identifying a largest sub-volume along the ray path and its corresponding maximum and minimum data values; checking if the opacity transfer function varies monotonically between the maximum and minimum scalar field values; if the function varies monotonically, calculating the amount of energy reflected by the sub-volume during its interaction with the ray; if the function does not vary monotonically, recursively identifying a smaller sub-volume along the ray path and its corresponding maximum and minimum data values; and checking if the opacity transfer function varies monotonically between the maximum and minimum scalar field values of the smaller sub-volume until a smaller sub-volume is identified in which the opacity transfer function varies monotonically or the smallest sub-volume is reached; and wherein the step of selectively casting a plurality of rays from a 2-D image plane includes subdividing the 2-D image plane into a plurality of sub-planes; for each of the plurality of sub-planes, casting four rays from the four corners of the sub-plane and estimating the pixel value of the reflected ray energy at each corner; calculating a maximum pixel value variation within the sub-plane; and recursively subdividing the sub-plane into multiple child sub-planes of smaller sizes by casting a ray from the center of the sub-plane until the maximum pixel value variation of the sub-plane is below a predefined imaging error threshold; and wherein the step of estimating the portion of the ray reflected by each sub-volume of the subset includes: estimating a maximum energy differential of the sub-volume; whereby said maximum energy differential measures the amount of approximation error introduced into reflected ray energy estimation when a sub-volume of varying transparency is replaced by another sub-volume of homogeneous transparency. comparing the maximum energy differential against a predefined energy error threshold; if the maximum energy differential is above the predefined energy error threshold, recursively selecting a smaller sub-volume along the ray path and estimating a new maximum energy differential of the smaller sub-volume; and if the maximum energy differential is below the predefined energy error threshold, calculating the amount of the ray energy reflected by the sub-volume using the illumination model. fragmenting a sampled 3-D dataset of a scalar field into a plurality of sub-volumes of different sizes, each sub-volume associated with a set of data value parameters characterizing the data value distribution of the scalar field within the sub-volume; defining an opacity transfer function that specifies a relationship between data values of the scalar field and opacity level; selectively casting a plurality of rays from a 2-D image plane towards the sampled dataset (106; 302), each ray having an initial ray energy and a cross-section; for each ray cast from a selected location on the 2-D image plane, adaptively selecting a subset of the plurality of sub-volumes for interaction with the ray based on the data value distribution of the scalar field within the sub-volumes that are located along the path of the ray; estimating the portion of the ray energy reflected by each sub-volume of the subset using an illumination model, the data values in the sub-volume and the opacity specified by the opacity transfer function for the data values in the sub-volume; and summing the reflected portions of the ray energy as a pixel value at the selected location on the 2-D image plane; and