Here we present a comparative computational study of the photoionization of furan from the ground and the two lowest-lying excited electronic states. The study aims to assess the quality of the computational methods currently employed for treating bound and continuum states in photoionization. For the ionization from the ground electronic state, we show that the Dyson orbital approach combined with an accurate solution of the continuum one particle wave functions in a multicenter B-spline basis, at the density functional theory (DFT) level, provides cross sections and asymmetry parameters in excellent agreement with experimental data. On the contrary, when the Dyson orbitals approach is combined with the Coulomb and orthogonalized Coulomb treatments of the continuum, the results are qualitatively different. In excited electronic states, three electronic structure methods, TDDFT, ADC(2), and CASSCF, have been used for the computation of the Dyson orbitals, while the continuum was treated at the B-spline/DFT level. We show that photoionization observables are sensitive probes of the nature of the excited states as well as of the quality of excited state wave functions. This paves the way for applications in more complex situations such as time resolved photoionization spectroscopy.