I'm using GDAL python to create a bunch of tiff files from a netcdf file.

def modelout_to_geotif(netcdfile, outdir, VarName):
    """ Given the path to the netcdf output,
    produces geotiffs for each band of the given variable name.
    
    Input::
        netcdfile: string path to netcdf file
        VarName: string variable name in netcdf file (e.g. 'GridPrec')
        
    Output::
        creates geotiffs for each 'band' (e.g. time step). Output is written
        to current working directory.
        
    
    """
    try:
        nci = gdal.Open('NETCDF:{0}:{1}'.format(netcdfile, VarName))
        ncd = Dataset(netcdfile)
    except:
        print "Could not open input file: {0}".format(netcdfile)
        return
    
    try:
        geotransform = nci.GetGeoTransform()
    except:
        print "Could not get geotransform.\n"
        return
    try:
        projection = nci.GetProjection()
    except:
        print "Could not get projection.\n"
        return
    
    numbands = nci.RasterCount
    times = [o2d(t) for t in ncd.variables['time']]
    print "Found {0} bands to process for {1}".format(numbands, VarName)
    for i in range(numbands):
        time = times[i]
        band = nci.GetRasterBand(i+1)
        raster = band.ReadAsArray()
        y,x = raster.shape
        
        output_file = "{0}_{1}.tif".format(VarName, time.strftime('%Y%m%d_%H'))
        subdir = os.path.join(outdir, VarName)
        if not os.path.exists(subdir):
            os.mkdir(subdir)
        output_file = os.path.join(subdir, output_file)
        if os.path.exists(output_file):
            continue
        # Create gtif
        driver = gdal.GetDriverByName("GTiff")
        #TODO: Make this more generic
        dst_ds = driver.Create(output_file, x, y, 1, gdal.GDT_CFloat32 )
        
        # top left x, w-e pixel resolution, rotation, top left y, rotation, n-s pixel resolution
        dst_ds.SetGeoTransform( geotransform )
  
        # set the reference info 
        dst_ds.SetProjection( projection )
        
        # write the band
        dst_ds.GetRasterBand(1).WriteArray(raster)

Now I would like to add those to a mosaic catalog in an ESRI File Geodatabase. I've managed to get it working using arcpy, but it's really slow, and I'd prefer to use GDAL. Is there a way? Below is my arcpy script:

import os
try:
    import arcpy
    from arcpy import env
except:
    raise ImportError("Unable to import arcpy - are you using ESRI Python?")

env = 'P:/Enki'


def add_to_mosaic(filegdb, raster_dir):    
    
    print env
    #Add Raster Dataset type Raster to FGDB Mosaic Dataset
    #Calculate Cell Size Ranges and Build Boundary
    #Build Overviews for Mosaic Dataset upon the 3rd level Raster Dataset pyramid
    #Apply TIFF file filter
    #Build Pyramids for the source datasets
    varnames = ['GridPrec']#,
                #'GridSnowOut',
                #'GridRad',
                #'GridTemp',
                #'GridSoilOut',
                #'GridRunoff']
    
    for VarName in varnames:
        gdbname = filegdb 
        mdname = os.path.join(gdbname, VarName)
        messages = 'Failed'
        sr =  arcpy.SpatialReference()    
#        sr.name = 'WGS_1984_UTM_Zone_32N'
        sr.factoryCode = 32632
        sr.create()
        print sr 
        try:
            noband = 1
            res = None
            #TODO: get from input data
            print 'trying to create %s' % VarName
            res = arcpy.CreateMosaicDataset_management(gdbname, mdname, sr, noband)

        except:
            print res.getMessages()
            
            
        
        
        rastype = "Raster Dataset"
        inpath = os.path.join(raster_dir, VarName)
        updatecs = "UPDATE_CELL_SIZES"
        updatebnd = "UPDATE_BOUNDARY"
        updateovr = "UPDATE_OVERVIEWS"
        maxlevel = "2"
        maxcs = "#"
        maxdim = "#"
        spatialref = "#"
        inputdatafilter = "*200701*.tif"
        subfolder = "NO_SUBFOLDERS"
        duplicate = "EXCLUDE_DUPLICATES"
        buildpy = "BUILD_PYRAMIDS"
        calcstats = "CALCULATE_STATISTICS"
        buildthumb = "NO_THUMBNAILS"
        comments = "Add Raster Datasets"
        forcesr = "#"
    
        res = arcpy.AddRastersToMosaicDataset_management(mdname,  rastype, inpath, updatecs, \
                                                   updatebnd, updateovr, maxlevel, maxcs, \
                                                   maxdim, spatialref, inputdatafilter, \
                                                   subfolder, duplicate, buildpy, calcstats, \
                                                   buildthumb, comments, forcesr)
        print res.getMessages()