The wave motion of a free field for general engineering can be simplified as a 1-D wave motion of an elastic layered half-space model，approximate solutions of which can be obtained by numerical methods.For such problems，the seismic waves are assumed to be vertical body waves propagating in the vertical direction, and site strata are regarded as nearly horizontal stratified structures.Even though there are many types of algorithms for seismic response analysis，all algorithms can be broadly classified into two main types.The first are numerical methods in the time domain，and the second are numerical methods in the frequency domain.However，numerical methods in the time domain，such as the finite difference method，finite element method，and boundary element method，are currently used.When these methods are used to calculate the free wave field in a layered half space，it is necessary to first discretize the computational region.The definite-solution problem of the continuous wave field is transformed to the problem of numerical computation of the discrete element nodes by methods of mathematical physics.The numerical formulas are usually expressed as a group of equations or explicit iteration schemes step-by-step in the time direction. However，the precision of approximate solutions computed by these numerical methods is affected by many factors，such as the mathematical algorithm，model range，mesh size，time step，and boundary condition.Inputting improper parameters will cause instability of the numerical algorithms，even causing no results to be obtained after a large amount of computation. Considering the generalized reflection transmission coefficient matrix method for synthetic seismograms，a new method is proposed，hich provides improvements for numerical methods in the time domain for solving problems of 1-D wave motion in an elastic layered half-space.When the method is used to compute a wave field in a layered half space，the element nodes are set at the wave impedance interfaces，which are called interfacial nodes.According to the wave motion principle of superposition，the wave field values between layers can be computed from the interface nodes，in which none of the nodes are set.Interfacial node values are made in accordance with the refracting and reflecting regulations at wave impedance interfaces and traveling time of waves between wave impedance interfaces，the expressions of which can be written as a group of time delay equations.Interfacial node values can be obtained after solving time delay equations.The wave field values in a layered half space can be obtained from one of the interface nodes at the free surface.Considering Huygens principle, the motion at interfacial nodes can be regarded as secondary sources or wavelet sources，when seismic waves pass through impedance interfaces.Therefore，the motion at interfacial nodes is called an interfacial wavelet，and the above method for determining the free wave field in a layered half space is called the interfacial wavelet method.The interfacial wavelet method is suitable for 1-D wave motion problems corresponding to the wave field of horizontal layered media caused by normal incidence.Two numerical results demonstrate that the proposed method has high accuracy and fast computing speed. In theory，the method can also be used to solve 2-D wave motion problems corresponding to the wave field of a horizontal layered media at oblique incidence.