The unprecedented resolution and sensitivity of ALMA have revolutionized our view of protostellar and protoplanetary disks, and our understanding of the star and planet formation process, and radiative transfer modeling has been a widely used tool helping to drive this revolution. The increasing power and availability of personal computing systems and supercomputers, however, as well as an increased focus on statistics and the development of suites of rigorous statistical tools suggest that we can improve the tools we use to study disks to maximize the scientific returns from our data. To that end I present pdspy, a publicly available code to do rigorous radiative transfer modeling of disks. pdspy links state-of-the-art Markov Chain Monte Carlo and Nested Sampling fitting with detailed dust and gas radiative transfer codes to rigorously model the structure of disks. It can handle a broad range of data, including continuum and spectral line millimeter observations, broadband spectra, and multi-wavelength imaging, and employs a flexible disk+envelope model with a large selection of tunable parameters. Here I will present the pdspy code, examples of science that can and have been done, and planned future capabilities in hopes that the community will find the code useful for their own science.