Skip to content

Most efficient metadata retrieval method #26

@joaomamede

Description

@joaomamede

I use pims_nd2 when my experiment files are broken and then repaired as it's the only library that works.
Basically I deconvolve the experiment, and I have a ome metadata xml generation function that retrieves several things of the metadata to write a .ome.tif.

It is quite slow to get the timestamp, x,y,z position in this script.
What I do is to loop through the reader at every time point and retrieve the metadata.

Is there anyway to retrieve timepoint specific metadata in a more efficient way?
Basically, what is slow is the "writeplanes" subfunction:

def pims_nd2meta2OMEXML(reader, project=False, time_offset=0, maxT=None,
        visit=0, order='TZCYX', pixeltype='uint16',
        verbose=False, **kwargs):
#     from apeer_ometiff_library import omexmlClass
    import aicsimageio.vendor.omexml as omexmlClass
    import re
    import sys,os

    #Missing TODO:
    #<Image>,  Name = "ImageName"
    #Instrument ID and Detector ID and Objective
    #     frames.parser._raw_metadata.image_metadata[b'SLxExperiment'][b'wsCameraName']

    # Objective settings with Refractive Index
    #Pixels,
    #Channel Color = RGB###, EmissionWavelength, Name of Channel.
    #Plane  ExposureTime, Position X, Y, Z (Z is easy as it's in nd2reader metadata)


    def writeplanes(pixel, v=0, reader=None, SizeT=1, SizeZ=1, SizeC=1, SizeV=1, order='TZCYX', timesteps=None, verbose=verbose, **kwargs):
        reader.default_coords['m'] = v
        if order == 'TZCYX':
            pixel.DimensionOrder = omexmlClass.DO_XYCZT
            counter = 0

            # timesteps = nd2meta['t_ms']/1000.
            #I don-t need this if maxT is smaller if doesn't reach the end of
            #the array anyway
            # if timesteps.shape[0] == SizeT*reader.sizes['z']:
            #     pass
            # else: timesteps = timesteps[:SizeT*reader.sizes['z']]

            for t in range(SizeT):
                for z in range(SizeZ):
                    for c in range(SizeC):
                        if verbose:
                            from IPython.display import clear_output
                            clear_output(wait = True)
                            print('Write PlaneTable: ', t, z, c),
                            sys.stdout.flush()

                        pixel.Plane(counter).TheT = t
                        pixel.Plane(counter).TheZ = z
                        pixel.Plane(counter).TheC = c

                        #check basically because of triggered acquisition the arrays shouldn't have the size of "channel"
                        #nd2reader spits ms /1000
                        #okay, this changes from triggered to non triggered I need to find a way to check if it was a triggered exp
                        #Fix to skip V's when they are not in the file

#                         print('tvzc',t,v,z,c)
                        pixel.Plane(counter).DeltaT = timesteps[t,z] + time_offset

#                         print(timesteps[v,t,z],pixel.Plane(counter).DeltaT)
#                         z_coords = np.array(
#                             nd2meta['z_coordinates']).reshape((reader.sizes['t'],reader.sizes['v'],reader.sizes['z']))
#                         if nd2meta['z_coordinates'] != None:
                        pixel.Plane(counter).PositionZ = reader[t].metadata['z_um'][z]
                        pixel.Plane(counter).PositionX = reader[t].metadata['x_um']
                        pixel.Plane(counter).PositionY = reader[t].metadata['y_um']
                        if verbose:       
                            print("y,x,z:",pixel.Plane(counter).PositionY, 
                                  pixel.Plane(counter).PositionX,
                                  pixel.Plane(counter).PositionZ)
                        counter = counter + 1
        return pixel

    #make a metadata var
    nd2meta = reader.metadata

    reader.iter_axes = 't'  # 't' is the default already
    reader.bundle_axes = 'zyx'  # when 'z' is available, this will be default
    reader.default_coords['m'] = visit
    #need to add the m component for correct time metadata.


    timesteps_array = np.zeros((reader.sizes['t'],reader.sizes['z']))
    for i in range(reader.sizes['t']):
            timesteps_array[i,:] = reader[i].metadata['t_ms']/1000

    # timesteps = np.zeros((reader.sizes['t']))
    # for i in range(reader.sizes['t']):
    #     timesteps[i] = reader[i].metadata['t_ms'][0]/1000

#     extra_dict = fetch_extra_metadata(reader)
#     Series = nd2meta['fields_of_view'][-1]+1
    scalex = reader.calibration
    scaley = scalex

    if not project:
#         scalez = round(nd2meta['z_coordinates'][1]-nd2meta['z_coordinates'][0],3)
#         scalez = frames.parser._raw_metadata.image_metadata[b'SLxExperiment'][b'ppNextLevelEx'][b''][b'ppNextLevelEx'][b''][b'uLoopPars'][b'dZStep']
        # scalez = reader[0].metadata['z_um'][1]-reader[0].metadata['z_um'][0]
        scalez = reader.calibrationZ
        z_coords = reader[0].metadata['z_um']

    # pixeltype = 'uint16'
    dimorder = order
# print(a)
    omexml = omexmlClass.OMEXML()
#     omexml.image_count = 1
#     omexml.image_count = reader.sizes['v']
    omexml.image(0).Name = os.path.basename(reader.filename)
    p = omexml.image(0).Pixels

    p.SizeX = reader.sizes['x']
    p.SizeY = reader.sizes['y']
    try:
        p.SizeC = reader.sizes['c']
    except: p.SizeC = 1
    try:
        SizeV = reader.sizes['m']
    except: SizeV = 1
#     SizeV = reader.sizes['v']
    if maxT == None:
        p.SizeT = reader.sizes['t']
        maxT = p.SizeT
    else: p.SizeT = maxT

    p.PhysicalSizeX = np.float(scalex)
    p.PhysicalSizeY = np.float(scaley)

    if not project:
        p.PhysicalSizeZ = np.float(scalez)
    p.PixelType = pixeltype
    try:
        p.channel_count = reader.sizes['c']
    except:
        p.channel_count = 1





    #I am using separate files for each visit point
    #, if you want one tiff with all visit points (possibly good for panels)
    #you will need to update this section
    for c in range(p.SizeC):
        p.Channel(c).Name = reader.metadata['plane_'+str(c)]['name']
        # try to automate by wavelenght one day, this basically sets the colour
        # to be automatically shown in Fiji/Image/NIS-elements
        clr = {'miRFP670':  65535 , 'mirfp670':  65535 ,'AF647': 65535,'a647': 65535,'Cy5': 65535,'640 nm': 65535,'pqbp1-AF647': 65535,
               'farRED-EM': 65535, 'mirfp67-': 65535,'Cy5 (Em)': 65535,
               'mruby3' : -16776961,'mRuby3' : -16776961,'mRuby' : -16776961,'RED-EM' : -16776961,'555 nm' : -16776961,'TRITC': -16776961, 'Cy3': -16776961,
                           'FITc': 16711935,   'fitc': 16711935,'GFP': 16711935,'FITC': 16711935,'GREEN-EM': 16711935, '470 nm': 16711935,'FITC (Em)': 16711935,'Igfp': 16711935, 'AF488': 16711935,'pre-paGFP': 16711935,'pre-paGFP': 16711935,'post-PAGFP':-16776961,'PAGFP':-16776961,
               'DAPI': 65535, 'Cgas-DY405': 65535, 'DAPI (Em)': 65535,'igfp': 16711935,}
        p.Channel(c).Color = clr[p.Channel(c).Name]
        p.Channel(c).ChannelEmissionWavelength = reader.metadata['plane_'+str(c)]['emission_nm']
        #Get this from metadata or ExtraData
        # p.Channel(c).EmissionWavelength =
        if pixeltype == 'unit8':
            #We never do averaging, update this if you do
            p.Channel(c).SamplesPerPixel = 1
        if pixeltype == 'unit16':
            p.Channel(c).SamplesPerPixel = 1


    if project:
        p.SizeZ = 1
        p.plane_count = 1 * maxT * p.SizeC #* SizeV
        p = writeplanes(p, v=visit, reader = reader, SizeT=maxT, SizeZ=1, SizeC=p.SizeC, SizeV=reader.sizes['m'], timesteps=timesteps_array,  order=dimorder)
    else:
        p.SizeZ = reader.sizes['z']
        p.plane_count = p.SizeZ * maxT * p.SizeC #* SizeV
        p = writeplanes(p, v=visit, reader = reader, SizeT=maxT, SizeZ=p.SizeZ, SizeC=p.SizeC, SizeV=reader.sizes['m'], timesteps=timesteps_array, order=dimorder)

    p.populate_TiffData()
    # omexml.structured_annotations.add_original_metadata(omexmlClass.OM_SAMPLES_PER_PIXEL, str(p.SizeC))
#     print(p.SizeT, p.SizeC, p.SizeZ)
#     print(omexml.to_xml())
    return omexml

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions