Combining far-infrared Herschel photometry from the PACS Evolutionary Probe (PEP) and Herschel Multi-tiered Extragalactic Survey (HerMES) guaranteed time programs with ancillary datasets in the GOODS-N, GOODS-S and COSMOS fields, it is possible to sample the 8-500 micron spectral energy distributions of galaxies with at least 7-10 bands. Extending to the UV, optical, and near- infrared, the number of bands increases up to 43. We reproduce the distribution of galaxies in a carefully selected 10 restframe color space, based on this rich data-set, using a superposition of multi-variate Gaussian modes. We use this model to classify galaxies and build median spectral energy distributions (SEDs) of each class, which are then fitted with a modified version of the MAGPHYS code that combines stellar light, emission from dust heated by stars and a possible warm dust contribution heated by an Active Galactic Nucleus (AGN). The color distribution of galaxies in each of the considered fields can be well described with the combination of 6-9 classes, spanning a large range of far- to near-IR luminosity ratios, as well as different strength of the AGN contribution to bolometric luminosities. The defined Gaussian grouping is used to identify rare or odd sources. The zoology of outliers includes Herschel-detected ellipticals, very blue z~1 Ly{alpha}-break galaxies, quiescent spirals, and torus-dominated AGN with star formation. Out of these groups and outliers, a new template library is assembled, consisting of 32 SEDs describing the intrinsic scatter in the restframe UV-to-submm colors of infrared galaxies. This library is tested against L(IR) estimates with and without Herschel data included, and compared to eight other popular methods often adopted in the literature. When implementing Herschel photometry, these approaches produce L(IR) values consistent with each other within a median absolute deviation of 10-20%, the scatter being dominated more by fine tuning of the codes, rather than by the choice of SED templates. Finally, the library is used to classify 24 micron detected sources in PEP GOODS fields on the basis of AGN content, L(60)/L(100) color and L(160)/L(1.6) luminosity ratio. AGN appear to be distributed in M*-SFR along with all other galaxies, regardless of the amount of infrared luminosity they are powering, with the tendency to lie on the high SFR side of the "main sequence". The incidence of warmer star-forming sources grows for objects with higher specific star formation rates, and they tend to populate the "off-sequence" region of the M*-SFR-z space.