From a system-theoretical point of view and for a given configuration of wells, there are only a limited number of degrees of freedom in the inputoutput dynamics of a reservoir system. This means that a large number of combinations of the state variables (pressure and saturation values) are not actually controllable and observable from the wells, and accordingly, they are not affecting the input-output behavior of the system. In an earlier publication, we therefore proposed a control-relevant upscaling methodology that uniformly coarsens the reservoir. Here, we present a control-relevant selective (i.e. non-uniform) coarsening (CRSC) method, in which the criterion for grid size adaptation is based on ranking the grid block contributions to the controllability and observability of the reservoir system. This multi-level CRSC method is attractive for use in iterative procedures such as computer-assisted flooding optimization for a given configuration of wells. In contrast to conventional flowbased coarsening techniques our method is independent of the specific flow rates or pressures imposed at the wells. Moreover the system-theoretical norms employed in our method provide tight upper bounds to the 'input-output energy' of the fine and coarse systems. These can be used as an a priori error-estimate of the performance of the coarse model. We applied our algorithm to two numerical examples and found that it can accurately reproduce results from the corresponding fine-scale simulations, while significantly speeding up the simulation.