We present 3D scattering waves imaging procedure based on the use of the velocity model in the time domain. The imaging is done by special summation of multicomponent multi shot - multi offset seismic data with special weights. The weights of this kind of summation we compute by tracing Gaussian beams from the imaging point towards the acquisition. This procedure gives the possibility to construct selective images corresponding to different directions of incident and scattering energy propagation that allows reconstruction of scattering indicatrices. To test this approach, we developed a realistic 3D geological model, containing typical diffraction objects, like faults, caves and fracture corridors. Synthetic seismic data for this model have been used as input data to verify the resolution and information content of the imaging procedure.