def load_dihedralfile(dihedralfile, zero_based_numbering=False):
    """
    Load definition of dihedral from a text file, i.e. Loading the file

    # i     j     k     j
      1     2     3     4
      2     3     4     5

    Will return dihedral_idxs = [(0,1,2,3), (1,2,3,4)]

    If a comment line #zero_based_numbering is found at the beginning of the file,
    or parameter zero_based_numbering == True, atom indices will be zero-based

    #zero_based_numbering
    # i     j     k     j
      1     2     3     4
      2     3     4     5

    Returns dihedral_idxs = [(1,2,3,4), (2,3,4,5)]

    If a fifth and sixth number are given in the line, they will be recognized as the lower
    and upper range limit of the dihedral angle, i.e. Reading the file

    # dihedral definition by atom indices starting from 1
    # i     j     k     j   (range_low)   (range_high)
      1     2     3     4     -120            120
      2     3     4     5      -90            150

    will generate two dihedrals, the first dihedral (0,1,2,3) have a range limit [-120, 120],
    the second dihedral (1,2,3,4) have the range limit [-90, 150], both ends inclusive.

    Parameters
    ----------
    dihedralfile: str
        filename that contains the dihedral angle definition
    zero_based_numbering: bool, default False
        Setting to true means atom indices in the file are zero-based

    Returns
    -------
    dihedral_idxs: list of list of 4 integers
        dihedrals are defined in 4 atom indices, e.g. [[i0,j0,k0,l0], [i1,j1,k1,l1]]
        The indices are zero-based.
    dihedral_ranges: list of list of 2 numbers
        dihedral ranges should be either empty (no limit), or two numbers [low, high],
        e.g. [[-120, 120], [-90, 150]]
        requirements: low >= -180, high <= 360, low < high
    """
    dihedral_idxs = []
    dihedral_ranges = []
    one_based_explicit = False # flag for user explicitly set one_based_numbering
    with open(dihedralfile) as infile:
        for line in infile:
            line = line.strip()
            if not line: continue
            if line[0] == '#':
                comment = line[1:].strip().lower()
                if comment == 'zero_based_numbering':
                    zero_based_numbering = True
                    if one_based_explicit is True:
                        raise ValueError("Can not specify both zero_based_numbering and one_based_indexing.")
                elif comment == 'one_based_numbering':
                    one_based_explicit = True
                    if zero_based_numbering is True:
                        raise ValueError("Can not specify both zero_based_numbering and one_based_indexing.")
            else:
                ls = line.split()
                if len(ls) == 4:
                    dihedral_idxs.append([int(i)-1 for i in ls])
                elif len(ls) == 6:
                    # insert default range [-180, 180] for missing values above
                    for _ in range(len(dihedral_idxs) - len(dihedral_ranges)):
                        dihedral_ranges.append([-180, 180])
                    dihedral_idxs.append([int(i)-1 for i in ls[:4]])
                    dihedral_ranges.append([int(v) for v in ls[4:]])
                else:
                    raise ValueError(f'Input line can not be recognized, should be 4 or 6 numbers\n{line}')
            # insert default range [-180, 180] for missing values
            if len(dihedral_ranges) > 0:
                for _ in range(len(dihedral_idxs) - len(dihedral_ranges)):
                    dihedral_ranges.append([-180, 180])
    # conver dihedral indices if zero based
    if zero_based_numbering:
        dihedral_idxs = [[i+1 for i in d] for d in dihedral_idxs]
    # check all dihedrals valid (>= 0)
    assert all(i >= 0 for d in dihedral_idxs for i in d), f'Dihedral indices {dihedral_idxs} error, all should >= 0'
    # check all ranges valid [-180, 180]
    assert all(low >= -180 and high <= 360 and low < high for low, high in dihedral_ranges), f'Dihedral ranges {dihedral_ranges} mistaken, range should be within [-180, 180]'
    # check dihedral_idxs and dihedral_ranges have same length
    if dihedral_ranges != []:
        assert len(dihedral_idxs) == len(dihedral_ranges), f'Dihedral ranges {dihedral_ranges} length not consistent with dihedral idxs {dihedral_idxs}'
    return dihedral_idxs, dihedral_ranges