A method of detecting subsurface features via a seismic survey, comprising: obtaining, by one or more processors (1210(a-n)) of a data processing system, seismic data detected by a plurality of seismic data acquisition units (30, 205) separated from each other by at least a threshold distance on a ground surface, the seismic data comprising image trace data for each of the plurality of seismic data acquisition units (30, 205); retrieving, by the data processing system from a binning function data structure in memory, a sampling interval and a parameter; configuring, by the data processing system, a bandlimited binning function with the sampling interval and the parameter; applying, by the data processing system, the bandlimited binning function to a plurality of image traces of the image trace data to generate a bandlimited angle gather value for a bin in an angle gathers grid (800, 900); generating, by the data processing system, an image based on the angle gathers grid (800, 900); and providing, by the data processing system, the image for display via a display device, the image indicating at least one of a subsurface lithological formation and a hydrocarbon deposit, characterized in that the amplitude of the amplitude of the 3D angle gather for the binning function is given by the formula: G(x,y;φ,γ,z) = Σ n I n (z) sinc[φ-φ n (z)] sinc[γ-γ n (z)], wherein G(x,y;φ,γ,z) represents an amplitude or the 3D angle gather at a location in the grid (x,y;φ,γ,z); Σ n represents a summation from n=1 to n=N number of image traces obtained from a plurality of seismic data acquisition units; I n represents an amplitude of an image trace at a fixed location (x,y) for a depth (z); sinc[φ-φ n (z)] represents the bandlimited binning function configured for φ, φ referring to the subsurface reflection azimuth angle, and sinc[γ-γ n (z)] represents the bandlimited binning function configured for γ, γ referring to the subsurface reflection angle. obtaining, by one or more processors (1210(a-n)) of a data processing system, seismic data detected by a plurality of seismic data acquisition units (30, 205) separated from each other by at least a threshold distance on a ground surface, the seismic data comprising image trace data for each of the plurality of seismic data acquisition units (30, 205); retrieving, by the data processing system from a binning function data structure in memory, a sampling interval and a parameter; configuring, by the data processing system, a bandlimited binning function with the sampling interval and the parameter; applying, by the data processing system, the bandlimited binning function to a plurality of image traces of the image trace data to generate a bandlimited angle gather value for a bin in an angle gathers grid (800, 900); generating, by the data processing system, an image based on the angle gathers grid (800, 900); and providing, by the data processing system, the image for display via a display device, the image indicating at least one of a subsurface lithological formation and a hydrocarbon deposit, characterized in that the amplitude of the amplitude of the 3D angle gather for the binning function is given by the formula: G(x,y;φ,γ,z) = Σ n I n (z) sinc[φ-φ n (z)] sinc[γ-γ n (z)], wherein G(x,y;φ,γ,z) represents an amplitude or the 3D angle gather at a location in the grid (x,y;φ,γ,z); Σ n represents a summation from n=1 to n=N number of image traces obtained from a plurality of seismic data acquisition units; I n represents an amplitude of an image trace at a fixed location (x,y) for a depth (z); sinc[φ-φ n (z)] represents the bandlimited binning function configured for φ, φ referring to the subsurface reflection azimuth angle, and sinc[γ-γ n (z)] represents the bandlimited binning function configured for γ, γ referring to the subsurface reflection angle.