title: Speed up the writing of hundreds of 3D numpy arrays to an hdf5 file
body: <p>the application I am working on convert a proprietary tiff file format (Nikon nd2 file) containing multiple images (fields of view), planes (Z) and fluorescence channels to numpy arrays that are then saved in an HDF5 file. Usually A tipical dataset has 50 field of view (fov) each of which has 5 channels and each channel has 40 z-planes). The overall file in around 6 Gb.</p>
  
  <p>This is the code I wrote:</p>
  
  <h1>Steps:</h1>
  
  <p>0) Import of all the required libraries</p>
  
  <pre><code>import nd2reader as nd2
  from matplotlib import pyplot as plt
  import numpy as np
  import h5py as h5
  import itertools
  import ast
  import glob as glob
  from joblib import Parallel, delayed
  import time
  </code></pre>
  
  <p>1) Function for running the conversion of the nd2 file. 
  The conversion to numpy arrays is done using nd2reader a python program and is quick. 
  To reduce the number of loops and use list comprehension I make a list of tuples each containing the channel and the fov 
  Example:
  [('DAPI', 0),
  ('DAPI', 1)]
  where DAPI is the channel and fov is the number.</p>
  
  <p>NB: THe experiment channel list is a file containing a dictionary than match the channel (key) with the gene of interest (value).</p>
  
  <pre><code>def ConvertND2File(ND2file):
  
      ChannelFileName=ND2file.replace('.nd2','ChannelsInfo.txt')
  
  
      # Read the file with the channels and raise an error if the file is missing
      try:
          ExperimentChannelList = ast.literal_eval(open(ChannelFileName).read())
      except IOError:
          print("The file:", ChannelFileName, "with the channels dictionary is missing")
          raise
  
      DataFileName=ND2file.replace('.nd2','.h5')
      with h5.File(DataFileName, 'w') as DataFile:
          ImgRef=nd2.Nd2(ND2file)
          Channels_Fields=itertools.product(ImgRef.channels,ImgRef.fields_of_view)
  
          # Create the empty array that will contain the 3D image
          ImgStack=np.empty([len(ImgRef.z_levels),ImgRef.height,ImgRef.width])
  
          # Use list comprehension to save the 3D arrays of the fov for each channel
          _=[SaveImg(DataFile,ImgRef,ExperimentChannelList,ImgStack,*x) for x in Channels_Fields]
  </code></pre>
  
  <p>2) Function to combine the images in a 3D array that is then written in an HDF5 file. I use h5py. I write each of the 3D numpy array on disk right after are generated. </p>
  
  <pre><code>def SaveImg(DataFile,ImgRef,ExperimentChannelList,ImgStack,*args):
          channel=args[0]
          fov=args[1]
          for idx,image in enumerate(ImgRef.select(channels=channel,z_levels=ImgRef.z_levels,fields_of_view=fov)):
              ImgStack[idx,:,:]=image
          gene=ExperimentChannelList[channel]
          ChannelGroup=DataFile.require_group(gene)
          FovDataSet=ChannelGroup.create_dataset(str(fov), data=ImgStack,dtype=np.float64,compression="gzip")
  </code></pre>
  
  <p>3) Body of the script and joblib call for parallel processing of all the files in a directory. </p>
  
  <pre><code>if __name__=='__main__':
  
      # Run the
      # Directory where ND2 file is stored (Ex. User/Data/)
      WorkingDirectory=input('Enter the directory with the files to process (ex. /User/):  ')
      #WorkingDirectory='/Users/simone/Box Sync/test/ND2conversion/'
      NumberOfProcesses=int(input('Enter the number of processes to use:  '))
      #NumberOfProcesses=2
      FileExt='nd2'
      # Iterator with the name of the files to process
      FilesIter=glob.iglob(WorkingDirectory+'*.'+FileExt)  
  
      now = time.time()
      Parallel(n_jobs=NumberOfProcesses,verbose=5)(delayed(ConvertND2File)(ND2file) for ND2file in FilesIter)
      print("Finished in", time.time()-now , "sec")
  </code></pre>
  
  <h1>Running time</h1>
  
  <p>Total time for the conversion of two file of 5.9 Gb<br>
  [Parallel(n_jobs=2)]: Done   1 out of   2 | elapsed:  7.4min remaining:  7.4min<br>
  [Parallel(n_jobs=2)]: Done   2 out of   2 | elapsed:  7.4min finished<br>
  Finished in 444.8717038631439 sec  </p>
  
  <h1>Question:</h1>
  
  <p>I was just wondering if there is a better way to handle the io to hdf5 files in order to speed up the conversion considering that if I want to scale up the process I will not be able to keep in memory all the 3D numpy arrays (fov) and then write them after each channel is processed.
  Thanks!</p>