def calc_system_ndof(natoms, nmolecs, nbonds, stop_com_tra, stop_com_rot):
    """
    Calculates the total / translational / rotational & internal /
    rotational / internal degrees of freedom of the system.

    Parameters
    ----------
    natoms : int
        Total number of atoms in the system
    nmolecs : int
        Total number of molecules in the system
    nbonds : int
        Total number of bonds in the system
    stop_com_tra : bool
        Was the center-of-mass translation removed during the simulation?
    stop_com_rot : bool
        Was the center-of-mass translation removed during the simulation?

    Returns
    -------
    ndof : dict
        Dictionary containing the degrees of freedom.
        Keys: ['tot', 'tra', 'rni', 'rot', 'int']
    """
    ndof_tot = 3 * natoms - nbonds
    if stop_com_tra:
        ndof_tot -= 3
    if stop_com_rot:
        ndof_tot -= 3
    ndof_tot_tra = 3 * nmolecs
    if stop_com_tra:
        ndof_tot -= 3
    ndof_tot_rni = ndof_tot - ndof_tot_tra
    ndof_tot_rot = 3 * nmolecs
    if stop_com_tra:
        ndof_tot -= 3
    ndof_tot_int = ndof_tot_rni - ndof_tot_rot
    ndof = {'tot': ndof_tot, 'tra': ndof_tot_tra, 'rni': ndof_tot_rni,
        'rot': ndof_tot_rot, 'int': ndof_tot_int}
    return ndof
