Imports a file in the Siemens .IMA format.
    :param filename: The filename of the file to import

def load_siemens_dicom(filename):
    """
    Imports a file in the Siemens .IMA format.
    :param filename: The filename of the file to import
    """
    # the .IMA format is a DICOM standard, unfortunately most of the information is contained inside a private and very
    # complicated header with its own data storage format, we have to get that information out along with the data
    # start by reading in the DICOM file completely
    dataset = pydicom.dicomio.read_file(filename)
    # now look through the tags (0029, 00xx) to work out which xx refers to the csa header
    # xx seems to start at 10 for Siemens
    xx = 0x0010
    header_index = 0
    while (0x0029, xx) in dataset:
        if dataset[0x0029, xx].value == "SIEMENS CSA HEADER":
            header_index = xx
        xx += 1
    # check that we have found the header
    if header_index == 0:
        raise KeyError("Could not find header index")
    # now we know which tag contains the CSA image header info: (0029, xx10)
    csa_header_bytes = dataset[0x0029, 0x0100 * header_index + 0x0010].value
    csa_header = read_csa_header(csa_header_bytes)
    #for key, value in csa_header.items():
    #    print("%s : %s" % (str(key), str(value)))
    # we can also get the series header info: (0029, xx20), but this seems to be mostly pretty boring

    # now we can work out the shape of the data (slices, rows, columns, fid_points)
    data_shape = (csa_header["SpectroscopyAcquisitionOut-of-planePhaseSteps"],
                  csa_header["Rows"],
                  csa_header["Columns"],
                  csa_header["DataPointColumns"],
                  )

    # now look through the tags (0029, 00xx) to work out which xx refers to the csa header
    # xx seems to start at 10 for Siemens
    xx = 0x0010
    data_index = 0
    while (0x7fe1, xx) in dataset:
        if dataset[0x7fe1, xx].value == "SIEMENS CSA NON-IMAGE":
            data_index = xx
        xx += 1
    # check that we have found the data
    if data_index == 0:
        raise KeyError("Could not find data index")
    # extract the actual data bytes
    csa_data_bytes = dataset[0x7fe1, 0x0100 * data_index + 0x0010].value
    # the data is stored as a list of 4 byte floats in (real, imaginary) pairs
    data_floats = struct.unpack("<%df" % (len(csa_data_bytes) / 4), csa_data_bytes)
    # form an iterator which will provide the numbers one at a time, then an iterator which calls that iterator twice
    # each cycle to give a complex pair
    data_iter = iter(data_floats)
    complex_iter = (complex(r, i) for r, i in zip(data_iter, data_iter))
    # give this iterator to numpy to build the data array
    complex_data = numpy.fromiter(complex_iter, "complex128", int(len(csa_data_bytes) / 8))
    # reshape the array to structure of rows, columns and slices
    complex_data = numpy.reshape(complex_data, data_shape).squeeze()

    return MRSData(complex_data,
                   csa_header["RealDwellTime"] * 1e-9,
                   csa_header["ImagingFrequency"],
                   te=csa_header["EchoTime"])