The purpose of this study was to implement the Monte Carlo method for clinical radiotherapy dose calculations. We used the EGS4/BEAM code to obtain the phase-space data for 6-20 MeV electron beams and 4, 6, and 15 MV photon beams for Varian Clinac 1800, 2100C, and 2300CD accelerators. A multiple-source model was used to reconstruct the phase-space data for both electron and photon beams, which retained the accuracy of the Monte Carlo beam data. The multiple-source model reduced the phase-space data storage requirement by a factor of 1000 and the accelerator simulation time by a factor of 10 or more. Agreement within 2% was achieved between the Monte Carlo calculations and measurements of the dose distributions in homogeneous and heterogeneous phantoms for various field sizes, source-surface distances, and beam modulations. The Monte Carlo calculated electron output factors were within 2% of the measured values for various treatment fields while the heterogeneity correction factors for various lung and bone phantoms were within 1% for photon beams and within 2% for electron beams. The EGS4/DOSXYZ Monte Carlo code was used for phantom and patient dose calculations. The results were compared to the dose distributions produced by a conventional treatment planning system and an intensity-modulated radiotherapy inverse-planning system. Significant differences (>5% in dose and >5 mm shift in isodose lines) were found between Monte Carlo calculations and the analytical calculations implemented in the commercial systems. Treatment sites showing the largest dose differences were for head and neck, lung, and breast cases.