A new approach based on a mixture of Gaussian and quadrilateral functions was developed to process bathymetric lidar waveforms. The approach was tested on two simulated data sets obtained from the existing Water-LIDAR (Wa-LID) waveform simulator. The first simulated data set corresponds to a sensor configuration modeled after a possible future satellite bathymetric lidar sensor that was previously studied. The second simulated data set corresponds to a lidar airborne configuration modeled using the HawkEye airborne lidar parameters. In the proposed approach, the lidar waveform is fitted into a combination of three functions, two Gaussians for both the water surface and water bottom contributions and a quadrilateral function to fit the water column contribution. The results show more accurate bathymetry estimates compared with the use of a triangular function to fit the column contribution or a simple peak detection method. For the satellite configuration, the bias is improved by 16.8 and 0.8 cm compared with the peak detection method and the use of a triangular function, respectively. For the airborne configuration, the bias is improved by 10.0 and 2.4 cm compared with the peak detection method and the use of a triangular function, respectively. The proposed waveform fitting using the quadrilateral function underestimates the bathymetry by -5.0 and - 6.1 cm for the simulated satellite and airborne data sets, respectively. The standard deviations of the bathymetry estimates are 6.0 and 8.2 cm, respectively. The obtained biases are inherent to overlaps between functions fitting the water surface, column, and bottom contributions.