Sorry for delayed reply, I was a little bit busy.
The general algorithm is the following.
1. Start with CIS or TDDFT calculation. Later you may use the obtained HF or DFT orbitals as starting approximation in CASSCF.
2. Look at the brightest state and find which orbitals contribute to it. In the attached example, the brightest state is #4, dominant contribution from MO 150->MO 154 transition, and three lower states are mostly MO 148,149,150->MO 151. So you should start with a small active space of MO 148, 149, 150, 151, 154. Just in case, I'd add MO 155, because is is quasidegenerate with MO 154 and is very similar to it.
3. Perform a state-averaged CASSCF with averaging over 5 states (ground state + 4 excited states).
4. Carefully examine what you obtain. In case of poor convergence, examine twice or even more. In this stage, you can vary the states to be averaged, active space, starting MOs.
5. Since CASSCF does not give transition intensities, you can perform TRANSITN run or better do a XMCQDPT calculation, which not only improves the energies, but gives transition moments as well. Thus, you have an absorption spectrum.
6. The absorption spectrum calculated at the ground-state geometry (optimized by RHF, MP2, DFT or state-specific CASSCF) should correspond to the experimental absorption spectrum. As for the emission, according to Kasha rule, emission occurs from the lowest excited state (with only few exceptions). So, to obtain an emission spectrum, you should always optimize S1 state. It is useful to turn on state tracking option in order to see, whether state flipping occurs during the optimization.
7. Calculate vertical spectrum at the optimized S1 geometry as indicated in p. 3, 4, 5.
>My biggest worry here is that I have no idea how to pick up CAS space for any inocuos molecule, e.g. something that does not have well known electronic transitions as in ample small molecule examples circulating in forum. Do I look at natural orbital occupation only in this case for values between ~1.9 and ~0.2? Looking at orbital shape in a large molecule almost doesn't matter (if it's not a bond breaking situation) as it doesn;t provide any information as those orbitals will hybridize over many atoms thus straying away from awell known conceptual orbitals.
Your're right, there is no simple prescription how to pick orbitals and states for averaging. See a good paper http://onlinelibrary.wiley.com/doi/10.1002/qua.23068/abstract
You should just try, and one day you get a good result.