def download_interim_for_gssha(main_directory,
                               start_datetime,
                               end_datetime,
                               leftlon=-180,
                               rightlon=180,
                               toplat=90,
                               bottomlat=-90,
                               precip_only=False):
    """
    Function to download ERA5 data for GSSHA

    .. note:: https://software.ecmwf.int/wiki/display/WEBAPI/Access+ECMWF+Public+Datasets

    Args:
        main_directory(:obj:`str`): Location of the output for the forecast data.
        start_datetime(:obj:`str`): Datetime for download start.
        end_datetime(:obj:`str`): Datetime for download end.
        leftlon(Optional[:obj:`float`]): Left bound for longitude. Default is -180.
        rightlon(Optional[:obj:`float`]): Right bound for longitude. Default is 180.
        toplat(Optional[:obj:`float`]): Top bound for latitude. Default is 90.
        bottomlat(Optional[:obj:`float`]): Bottom bound for latitude. Default is -90.
        precip_only(Optional[bool]): If True, will only download precipitation.

    Example::

        from gsshapy.grid.era_to_gssha import download_era_interim_for_gssha

        era_interim_folder = '/era_interim'
        leftlon = -95
        rightlon = -75
        toplat = 35
        bottomlat = 30
        download_era_interim_for_gssha(era5_folder, leftlon, rightlon, toplat, bottomlat)

    """
    # parameters: https://software.ecmwf.int/wiki/display/CKB/Details+of+ERA-Interim+parameters

    # import here to make sure it is not required to run
    from ecmwfapi import ECMWFDataServer
    server = ECMWFDataServer()

    try:
        mkdir(main_directory)
    except OSError:
        pass

    download_area = "{toplat}/{leftlon}/{bottomlat}/{rightlon}".format(toplat=toplat,
                                                               leftlon=leftlon,
                                                               bottomlat=bottomlat,
                                                               rightlon=rightlon)
    download_datetime = start_datetime
    interim_request = {
        'dataset': "interim",
        #  'oper' specifies the high resolution daily data, as opposed to monthly means, wave, eda edmm, etc.
        'stream': "oper",
        #  Surface level, as opposed to pressure level (pl) or model level (ml)
        'levtype': "sfc",
        # The spatial resolution in ERA interim is 80 km globally on a Gaussian grid.
        # Here we us lat/long with 0.75 degrees, which is approximately the equivalent of 80km.
        'grid': "0.5/0.5",
        'area': download_area,
        'format': 'netcdf',
    }
    while download_datetime <= end_datetime:
        interim_request['date'] = download_datetime.strftime("%Y-%m-%d")
        if not precip_only:
            download_file = path.join(main_directory, "erai_gssha_{0}_an.nc".format(download_datetime.strftime("%Y%m%d")))
            if not path.exists(download_file):
                #  We want instantaneous parameters, which are archived as type Analysis ('an') as opposed to forecast (fc)
                interim_request['type'] = "an"
                # For parameter codes see the ECMWF parameter database at http://apps.ecmwf.int/codes/grib/param-db
                interim_request['param'] = "2t/2d/sp/10u/10v/tcc"
                # step 0 is analysis, 3-12 is forecast
                interim_request['step'] = "0"
                # ERA Interim provides 6-hourly analysis
                interim_request['time'] = "00/06/12/18"
                interim_request['target'] = download_file
                server.retrieve(interim_request)

            download_file = path.join(main_directory, "erai_gssha_{0}_1_fc.nc".format(download_datetime.strftime("%Y%m%d")))
            if not path.exists(download_file):
                interim_request['type'] = "fc"
                interim_request['param'] = "2t/2d/sp/10u/10v/tcc"
                interim_request['step'] = "3"
                interim_request['time'] = "00/06/12/18"
                interim_request['target'] = download_file
                server.retrieve(interim_request)

        download_file = path.join(main_directory, "erai_gssha_{0}_fc.nc".format(download_datetime.strftime("%Y%m%d")))
        if not path.exists(download_file):
            interim_request['type'] = "fc"
            interim_request['param'] = "tp/ssrd"
            interim_request['step'] = "3/6/9/12"
            interim_request['time'] = "00/12"
            interim_request['target'] = download_file
            server.retrieve(interim_request)
            # TODO: READ FILE AND MODIFY VALUES SO IT IS NOT INCREMENTAL
            # https://software.ecmwf.int/wiki/pages/viewpage.action?pageId=56658233
            # You need  total precipitation for every 6 hours.
            # Daily total precipitation (tp) is only available with a forecast base time 00:00 and 12:00,
            # so to get tp for every 6 hours you will need to extract (and for the second and fourth period calculate):
            # tp(00-06) = (time 00, step 6)
            # tp(06-12) = (time 00, step 12) minus (time 00, step 6)
            # tp(12-18) = (time 12, step 6)
            # tp(18-24) = (time 12, step 12) minus (time 12, step 6)
            # (Note the units for total precipitation is meters.)
            tmp_download_file = download_file + '_tmp'
            with xr.open_dataset(download_file) as xd:
                diff_xd = xd.diff('time')
                xd.tp[1:4] = diff_xd.tp[:3]
                xd.tp[5:] = diff_xd.tp[4:]
                xd.ssrd[1:4] = diff_xd.ssrd[:3]
                xd.ssrd[5:] = diff_xd.ssrd[4:]
                xd.to_netcdf(tmp_download_file)
            remove(download_file)
            rename(tmp_download_file, download_file)

        download_file = path.join(main_directory, "erai_gssha_{0}_0_fc.nc".format(download_datetime.strftime("%Y%m%d")))
        if download_datetime <= start_datetime and not path.exists(download_file):
            loc_download_date = (download_datetime-timedelta(1)).strftime("%Y-%m-%d")
            interim_request['type'] = "fc"
            interim_request['param'] = "tp/ssrd"
            interim_request['step'] = "9/12"
            interim_request['time'] = "12"
            interim_request['target'] = download_file
            interim_request['date'] = loc_download_date
            server.retrieve(interim_request)
            # convert to incremental (see above)
            tmp_download_file = download_file + '_tmp'
            with xr.open_dataset(download_file) as xd:
                inc_xd = xd.diff('time')
                inc_xd.to_netcdf(tmp_download_file)
            remove(download_file)
            rename(tmp_download_file, download_file)
        download_datetime += timedelta(1)