Terahertz Time Domain Spectroscopy (THz-TDS) systems have emerged as mature technologies with significant potential across various research fields and industries. However, the lack of standardized methods for signal and noise estimation and reduction hinders its full potential. This paper introduces a methodology to significantly reduce noise in THz-TDS time traces, providing a reliable and less biased estimation of the signal. The method results in an improved signal-to-noise ratio, enabling the utilization of the full dynamic range of such setups. Additionally, we investigate the estimation of the covariance matrix to quantify the uncertainties associated with the signal estimator. This matrix is essential for extracting accurate material parameters by normalizing the error function in the fitting process. Our approach addresses practical scenarios where the number of repeated measurements is limited compared to the sampling time axis length. We envision this work as the initial step toward standardizing THz-TDS data processing. We believe it will foster collaboration between the THz and signal processing communities, leading to the development of more sophisticated methods to tackle new challenges introduced by novel setups based on optoelectronic devices and dual-comb spectroscopy.