We report on a first-principle theoretical investigation of the optical absorption and emission spectra of poly(3-hexylthiophene) (P3HT) aggregates by means of a multiscale all-atom hybrid approach, which combines: (i) molecular dynamics simulations of (intrachain) conformational and (interchain) positional disorder, (ii) quantum-chemical calculations of intrachain excited states and excitonic interactions, and (iii) solving of a Frenkel–Holstein model that includes coupling of the electronic excitations to a dominant high-frequency molecular vibration. The modeling work points to the co-existence of ordered nano-aggregates into disordered domains, with the former dominating the emission spectrum while both regions contribute to the absorption spectrum. The measured photoluminescence line shape and the Stokes shift are both reproduced by the theory for aggregates comprising polymer chain with lengths ∼35–40 repeating units.