This study aims to build a 3D velocity model for investigating and evaluating underground contaminants flowing through groundwater pathways. To accomplish this goal, we performed seismic exploration at a test site with fractured rock layers within a depth of 100 m. The acquired seismic data was then used to construct a three-dimensional (3D) P-wave velocity model of the test site. Although comprehensive geophysical exploration may be useful for investigating the near-surface structure, this paper focuses only on constructing a P-wave velocity model using the seismic method. The primary information of the model consists of extracted velocities from a two-dimensional surface, borehole first-arrival traveltime tomography results, and full waveform sonic log data. Since the test site has the spatial restriction of the survey line, we intended to improve the geological structure analysis results using various quantitative and qualitative analysis methods. First, to increase the reliability of velocity information, we performed a traveltime analysis on the zero vertical interval (ZVI) gather from the borehole seismic data. Then, we identified the qualitative information of the fracture zone's location by analyzing the amplitude variation of the ZVI gather. Moreover, we extracted structural information using the common reflection point gather from the borehole seismic data to supplement the obtained lithological information. However, the information for constructing a 3D P-wave velocity model was still insufficient due to the spatial constraints of the survey line and the limited depth of the borehole seismic survey. We filled this gap using the radial basis function interpolation method. We could verify the completed 3D P-wave velocity model by comparing it with core log data. Overall, the integrated interpretation of the final velocity model and the analysis results could provide a probable pathway that indicates hydraulic connectivity.