We present the results of an extensive lattice calculation of the renormalization constants of bilinear and four-quark operators for the non-perturbatively O(a)improved Wilson action. The results are obtained in the quenched approximation at four values of the lattice coupling by using the non-perturbative RI/MOM renormalization method. Several sources of systematic uncertainties, including discretization errors and final volume effects, are examined. The contribution of the Goldstone pole, which in some cases may affect the extrapolation of the renormalization constants to the chiral limit, is non-perturbatively subtracted. The scale independent renormalization constants of bilinear quark operators have been also computed by using the lattice chiral Ward identities approach and compared with those obtained with the RI-MOM method. For those renormalization constants the non-perturbative estimates of which have been already presented in the literature we find an agreement which is typically at the level of 1%.PACS numbers: 11.15.Ha, 12.38.Gc, 11.10.Gh.1 See also ref.[7] for a related study in the context of the lattice non-linear σ-model.• Finite volume effects have been examined by comparing the results of two independent simulations performed, at the same value of the lattice coupling (β = 6.0), on different lattice volumes.• Goldstone pole contributions: power suppressed contributions coming from the Goldstone pole, which may affect the extrapolation to the chiral limit of the RCs of the pseudoscalar density and of the four-fermion operators coupled to the Goldstone boson, have been non-perturbatively subtracted.The scale independent RCs of the vector and axial-vector currents, Z V and Z A , and the ratio Z P /Z S , have also been determined in this study by using the lattice chiral WI approach, and the results are compared with those obtained with the RI/MOM method. We also compare our determinations of these constants with those obtained by the AL-PHA [16,17] and LANL [18] Collaborations by using the WI method, and our determination of Z P with the one obtained by ALPHA [19] within the SF approach. We find an agreement which is typically at the level of 1%. Since the systematics involved in these approaches are different, these comparisons provide additional confidence on the high level of accuracy reached in the implementation of the RI/MOM non-perturbative method.Our final results for the RCs of bilinear quark operators and four-fermion operators are collected in tables 2-7. Preliminary results of the present study have been presented at the Lattice 2002 conference [5,20].The plan of this paper is as follows. In sec. 2 we give the details of the numerical simulation and briefly review the non-perturbative RI/MOM method used to determine lattice RCs. The various sources of systematic errors involved in the calculation are discussed in sec. 3, where we also provide estimates of the corresponding uncertainties. The results for the RCs of bilinear quark operators obtained with the RI/MOM method are summarized i...