Source code for xmipp3.convert.convert

# **************************************************************************
# *
# * Authors:     J.M. De la Rosa Trevin (jmdelarosa@cnb.csic.es)
# *              Laura del Cano (ldelcano@cnb.csic.es)
# *
# * Unidad de  Bioinformatica of Centro Nacional de Biotecnologia , CSIC
# *
# * This program is free software; you can redistribute it and/or modify
# * it under the terms of the GNU General Public License as published by
# * the Free Software Foundation; either version 2 of the License, or
# * (at your option) any later version.
# *
# * This program is distributed in the hope that it will be useful,
# * but WITHOUT ANY WARRANTY; without even the implied warranty of
# * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# * GNU General Public License for more details.
# *
# * You should have received a copy of the GNU General Public License
# * along with this program; if not, write to the Free Software
# * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
# * 02111-1307  USA
# *
# *  All comments concerning this program package may be sent to the
# *  e-mail address 'scipion@cnb.csic.es'
# *
# **************************************************************************
"""
This module contains converter functions that will serve to:
1. Write from base classes to Xmipp specific files
2. Read from Xmipp files to base classes
"""

import os
from os.path import join, exists
from collections import OrderedDict
try:
    from itertools import izip
except ImportError:
    izip = zip
import numpy as np

from pyworkflow.utils import replaceBaseExt
from pyworkflow.object import ObjectWrap, String, Float, Integer, Object
from pwem.constants import (NO_INDEX, ALIGN_NONE, ALIGN_PROJ, ALIGN_2D,
                            ALIGN_3D)
from pwem.objects import (Angles, Coordinate, Micrograph, Volume, Particle,
                          MovieParticle, CTFModel, Acquisition, SetOfParticles,
                          Class3D, SetOfVolumes, Transform)
from pwem.emlib.image import ImageHandler
import pwem.emlib.metadata as md

from xmipp3.base import getLabelPythonType, iterMdRows

from pwem import emlib

[docs]def prefixAttribute(attribute): return '_xmipp_%s' % attribute
if not getattr(emlib, "GHOST_ACTIVATED", False): """ Some of MDL may not exist when Ghost is activated """ # This dictionary will be used to map # between CTFModel properties and Xmipp labels ACQUISITION_DICT = OrderedDict([ ("_amplitudeContrast", emlib.MDL_CTF_Q0), ("_sphericalAberration", emlib.MDL_CTF_CS), ("_voltage", emlib.MDL_CTF_VOLTAGE) ]) COOR_DICT = OrderedDict([ ("_x", emlib.MDL_XCOOR), ("_y", emlib.MDL_YCOOR) ]) COOR_EXTRA_LABELS = [ # Additional autopicking-related metadata md.RLN_PARTICLE_AUTOPICK_FOM, md.RLN_PARTICLE_CLASS, md.RLN_ORIENT_PSI, emlib.MDL_GOOD_REGION_SCORE ] CTF_DICT = OrderedDict([ ("_defocusU", emlib.MDL_CTF_DEFOCUSU), ("_defocusV", emlib.MDL_CTF_DEFOCUSV), ("_defocusAngle", emlib.MDL_CTF_DEFOCUS_ANGLE), ("_resolution", emlib.MDL_CTF_CRIT_MAXFREQ), ("_fitQuality", emlib.MDL_CTF_CRIT_FITTINGSCORE) ]) # TODO: remove next dictionary when all # cTFmodel has resolution and fitQuality CTF_DICT_NORESOLUTION = OrderedDict([ ("_defocusU", emlib.MDL_CTF_DEFOCUSU), ("_defocusV", emlib.MDL_CTF_DEFOCUSV), ("_defocusAngle", emlib.MDL_CTF_DEFOCUS_ANGLE) ]) CTF_PSD_DICT = OrderedDict([ ("_psdFile", emlib.MDL_PSD), (prefixAttribute("enhanced_psd"), emlib.MDL_PSD_ENHANCED), (prefixAttribute("ctfmodel_quadrant"), emlib.MDL_IMAGE1), (prefixAttribute("ctfmodel_halfplane"), emlib.MDL_IMAGE1) ]) CTF_EXTRA_LABELS = [ emlib.MDL_CTF_CA, emlib.MDL_CTF_ENERGY_LOSS, emlib.MDL_CTF_LENS_STABILITY, emlib.MDL_CTF_CONVERGENCE_CONE, emlib.MDL_CTF_LONGITUDINAL_DISPLACEMENT, emlib.MDL_CTF_TRANSVERSAL_DISPLACEMENT, emlib.MDL_CTF_K, emlib.MDL_CTF_BG_GAUSSIAN_K, emlib.MDL_CTF_BG_GAUSSIAN_SIGMAU, emlib.MDL_CTF_BG_GAUSSIAN_SIGMAV, emlib.MDL_CTF_BG_GAUSSIAN_CU, emlib.MDL_CTF_BG_GAUSSIAN_CV, emlib.MDL_CTF_BG_SQRT_K, emlib.MDL_CTF_BG_SQRT_U, emlib.MDL_CTF_BG_SQRT_V, emlib.MDL_CTF_BG_SQRT_ANGLE, emlib.MDL_CTF_BG_BASELINE, emlib.MDL_CTF_BG_GAUSSIAN2_K, emlib.MDL_CTF_BG_GAUSSIAN2_SIGMAU, emlib.MDL_CTF_BG_GAUSSIAN2_SIGMAV, emlib.MDL_CTF_BG_GAUSSIAN2_CU, emlib.MDL_CTF_BG_GAUSSIAN2_CV, emlib.MDL_CTF_BG_GAUSSIAN2_ANGLE, emlib.MDL_CTF_CRIT_FITTINGCORR13, emlib.MDL_CTF_CRIT_ICENESS, emlib.MDL_CTF_VPP_RADIUS, emlib.MDL_CTF_DOWNSAMPLE_PERFORMED, emlib.MDL_CTF_CRIT_PSDVARIANCE, emlib.MDL_CTF_CRIT_PSDPCA1VARIANCE, emlib.MDL_CTF_CRIT_PSDPCARUNSTEST, emlib.MDL_CTF_CRIT_FIRSTZEROAVG, emlib.MDL_CTF_CRIT_DAMPING, emlib.MDL_CTF_CRIT_FIRSTZERORATIO, emlib.MDL_CTF_CRIT_PSDCORRELATION90, emlib.MDL_CTF_CRIT_PSDRADIALINTEGRAL, emlib.MDL_CTF_CRIT_NORMALITY, # In xmipp the ctf also contains acquisition information emlib.MDL_CTF_Q0, emlib.MDL_CTF_CS, emlib.MDL_CTF_VOLTAGE, emlib.MDL_CTF_SAMPLING_RATE ] # TODO: remove next dictionary when all # cTFmodel has resolution and fitquality CTF_EXTRA_LABELS_PLUS_RESOLUTION = [ emlib.MDL_CTF_CA, emlib.MDL_CTF_ENERGY_LOSS, emlib.MDL_CTF_LENS_STABILITY, emlib.MDL_CTF_CONVERGENCE_CONE, emlib.MDL_CTF_LONGITUDINAL_DISPLACEMENT, emlib.MDL_CTF_TRANSVERSAL_DISPLACEMENT, emlib.MDL_CTF_K, emlib.MDL_CTF_BG_GAUSSIAN_K, emlib.MDL_CTF_BG_GAUSSIAN_SIGMAU, emlib.MDL_CTF_BG_GAUSSIAN_SIGMAV, emlib.MDL_CTF_BG_GAUSSIAN_CU, emlib.MDL_CTF_BG_GAUSSIAN_CV, emlib.MDL_CTF_BG_SQRT_K, emlib.MDL_CTF_BG_SQRT_U, emlib.MDL_CTF_BG_SQRT_V, emlib.MDL_CTF_BG_SQRT_ANGLE, emlib.MDL_CTF_BG_BASELINE, emlib.MDL_CTF_BG_GAUSSIAN2_K, emlib.MDL_CTF_BG_GAUSSIAN2_SIGMAU, emlib.MDL_CTF_BG_GAUSSIAN2_SIGMAV, emlib.MDL_CTF_BG_GAUSSIAN2_CU, emlib.MDL_CTF_BG_GAUSSIAN2_CV, emlib.MDL_CTF_BG_GAUSSIAN2_ANGLE, emlib.MDL_CTF_CRIT_MAXFREQ, # ### emlib.MDL_CTF_CRIT_FITTINGSCORE, # ### emlib.MDL_CTF_CRIT_FITTINGCORR13, emlib.MDL_CTF_CRIT_ICENESS, emlib.MDL_CTF_DOWNSAMPLE_PERFORMED, emlib.MDL_CTF_CRIT_PSDVARIANCE, emlib.MDL_CTF_CRIT_PSDPCA1VARIANCE, emlib.MDL_CTF_CRIT_PSDPCARUNSTEST, emlib.MDL_CTF_CRIT_FIRSTZEROAVG, emlib.MDL_CTF_CRIT_DAMPING, emlib.MDL_CTF_CRIT_FIRSTZERORATIO, emlib.MDL_CTF_CRIT_PSDCORRELATION90, emlib.MDL_CTF_CRIT_PSDRADIALINTEGRAL, emlib.MDL_CTF_CRIT_NORMALITY, # In xmipp the ctf also contains acquisition information emlib.MDL_CTF_Q0, emlib.MDL_CTF_CS, emlib.MDL_CTF_VOLTAGE, emlib.MDL_CTF_SAMPLING_RATE, emlib.MDL_CTF_VPP_RADIUS, ] # Some extra labels to take into account the zscore IMAGE_EXTRA_LABELS = [ emlib.MDL_ZSCORE, emlib.MDL_ZSCORE_HISTOGRAM, emlib.MDL_ZSCORE_RESMEAN, emlib.MDL_ZSCORE_RESVAR, emlib.MDL_ZSCORE_RESCOV, emlib.MDL_ZSCORE_SHAPE1, emlib.MDL_ZSCORE_SHAPE2, emlib.MDL_ZSCORE_SNR1, emlib.MDL_ZSCORE_SNR2, emlib.MDL_CUMULATIVE_SSNR, emlib.MDL_PARTICLE_ID, emlib.MDL_FRAME_ID, emlib.MDL_SCORE_BY_VAR, emlib.MDL_SCORE_BY_GINI, emlib.MDL_ZSCORE_DEEPLEARNING1 ] ANGLES_DICT = OrderedDict([ ("_angleY", emlib.MDL_ANGLE_Y), ("_angleY2", emlib.MDL_ANGLE_Y2), ("_angleTilt", emlib.MDL_ANGLE_TILT) ]) ALIGNMENT_DICT = OrderedDict([ (prefixAttribute("shiftX"), emlib.MDL_SHIFT_X), (prefixAttribute("shiftY"), emlib.MDL_SHIFT_Y), (prefixAttribute("shiftZ"), emlib.MDL_SHIFT_Z), (prefixAttribute("flip"), emlib.MDL_FLIP), (prefixAttribute("anglePsi"), emlib.MDL_ANGLE_PSI), (prefixAttribute("angleRot"), emlib.MDL_ANGLE_ROT), (prefixAttribute("angleTilt"), emlib.MDL_ANGLE_TILT), ])
[docs]def objectToRow(obj, row, attrDict, extraLabels=[]): """ This function will convert an EMObject into a Row. Params: obj: the EMObject instance (input) row: the Row instance (output) -see emlib.metadata.utils.Row()- attrDict: dictionary with the map between obj attributes(keys) and row MDLabels in Xmipp (values). extraLabels: a list with extra labels that could be included as _xmipp_labelName """ if obj.isEnabled(): enabled = 1 else: enabled = -1 row.setValue(emlib.MDL_ENABLED, enabled) for attr, label in attrDict.items(): if hasattr(obj, attr): valueType = getLabelPythonType(label) value = getattr(obj, attr).get() try: row.setValue(label, valueType(value)) except Exception as e: print(e) print("Problems found converting metadata: ") print("Label id = %s" % label) print("Attribute = %s" % attr) print("Value = %s" % value) print("Value type = %s" % valueType) raise e row.setValue(label, valueType(getattr(obj, attr).get())) attrLabels = attrDict.values() for label in extraLabels: attrName = prefixAttribute(emlib.label2Str(label)) if label not in attrLabels and hasattr(obj, attrName): value = obj.getAttributeValue(attrName) row.setValue(label, value)
[docs]def rowToObject(row, obj, attrDict, extraLabels=[]): """ This function will convert from a Row to an EMObject. Params: row: the Row instance (input) -see emlib.metadata.utils.Row()- obj: the EMObject instance (output) attrDict: dictionary with the map between obj attributes(keys) and row MDLabels in Xmipp (values). extraLabels: a list with extra labels that could be included as _xmipp_labelName """ obj.setEnabled(row.getValue(emlib.MDL_ENABLED, 1) > 0) for attr, label in attrDict.items(): value = row.getValue(label) if not hasattr(obj, attr): setattr(obj, attr, ObjectWrap(value)) else: getattr(obj, attr).set(value) attrLabels = attrDict.values() for label in extraLabels: if label not in attrLabels and row.hasLabel(label): labelStr = emlib.label2Str(label) setattr(obj, prefixAttribute(labelStr), row.getValueAsObject(label))
[docs]def setXmippAttributes(obj, objRow, *labels, **kwargs): """ Set the labels attribute(s) to obj from a objRow (those not basic ones). More than one label can be set and default can also be a list of same number of entries. The new attribute will be named _xmipp_LabelName and the datatype will be set correctly. If obj doesn't contain a label and a default value is in kwargs, it's set in its corresponding Scipion datatype. """ defaults = kwargs.get('default', None) for idx, label in enumerate(labels): if objRow.containsLabel(label): value = objRow.getValueAsObject(label) else: if isinstance(defaults, list): value = defaults[idx] else: value = defaults value = getScipionObj(value) if value is not None: setXmippAttribute(obj, label, value)
[docs]def getScipionObj(value): if isinstance(value, Object): return value elif isinstance(value, int): return Integer(value) elif isinstance(value, float): return Float(value) elif isinstance(value, str): return String(value) else: return None
[docs]def setXmippAttribute(obj, label, value): """ Sets an attribute of an object prefixing it with xmipp""" setattr(obj, prefixAttribute(emlib.label2Str(label)), value)
[docs]def getXmippAttribute(obj, label, default=None): """ Sets an attribute of an object prefixing it with xmipp""" return getattr(obj, prefixAttribute(emlib.label2Str(label)), default)
[docs]def rowFromMd(mdIn, objId): row = md.Row() row.readFromMd(mdIn, objId) return row
def _containsAll(row, labels): """ Check if the labels (values) in labelsDict are present in the row. """ values = labels.values() if isinstance(labels, dict) else labels return all(row.containsLabel(l) for l in values) def _containsAny(row, labels): """ Check if the labels (values) in labelsDict are present in the row. """ values = labels.values() if isinstance(labels, dict) else labels return any(row.containsLabel(l) for l in values) def _rowToObjectFunc(obj): """ From a given object, return the function rowTo{OBJECT_CLASS_NAME}. """ return globals()['rowTo' + obj.getClassName()] def _readSetFunc(obj): """ From a given object, return the function read{OBJECT_CLASS_NAME}. """ return globals()['read' + obj.getClassName()]
[docs]def locationToXmipp(index, filename): """ Convert an index and filename location to a string with @ as expected in Xmipp. """ return ImageHandler.locationToXmipp((index, filename))
[docs]def fixVolumeFileName(image): """ This method will add :mrc to .mrc volumes because for mrc format is not possible to distinguish between 3D volumes and 2D stacks. """ return ImageHandler.fixXmippVolumeFileName(image)
[docs]def getMovieFileName(movie): """ Add the :mrcs or :ems extensions to movie files to be recognized by Xmipp as proper stack files. """ fn = movie.getFileName() if fn.endswith('.mrc'): fn += ':mrcs' elif fn.endswith('.em'): fn += ':ems' return fn
[docs]def getImageLocation(image): return ImageHandler.locationToXmipp(image)
[docs]def xmippToLocation(xmippFilename): """ Return a location (index, filename) given a Xmipp filename with the index@filename structure. """ if '@' in xmippFilename: return emlib.FileName(xmippFilename).decompose() else: return NO_INDEX, str(xmippFilename)
[docs]def setObjId(obj, mdRow, label=emlib.MDL_ITEM_ID): if mdRow.containsLabel(label): obj.setObjId(mdRow.getValue(label)) else: obj.setObjId(None)
[docs]def setRowId(mdRow, obj, label=emlib.MDL_ITEM_ID): mdRow.setValue(label, int(obj.getObjId()))
[docs]def micrographToCTFParam(mic, ctfparam): """ This function is used to convert a Micrograph object to the .ctfparam metadata needed by some Xmipp programs. If the micrograph already comes from xmipp, the ctfparam will be returned, if not, the new file. """ ctf = mic.getCTF() if hasattr(ctf, '_xmippMd'): return ctf._xmippMd.get() mdCtf = emlib.MetaData() mdCtf.setColumnFormat(False) row = md.Row() ctfModelToRow(ctf, row) acquisitionToRow(mic.getAcquisition(), row) row.writeToMd(mdCtf, mdCtf.addObject()) mdCtf.write(ctfparam) return ctfparam
[docs]def imageToRow(img, imgRow, imgLabel, **kwargs): # Provide a hook to be used if something is needed to be # done for special cases before converting image to row preprocessImageRow = kwargs.get('preprocessImageRow', None) if preprocessImageRow: preprocessImageRow(img, imgRow) setRowId(imgRow, img) # Set the id in the metadata as MDL_ITEM_ID index, filename = img.getLocation() fn = locationToXmipp(index, filename) imgRow.setValue(imgLabel, fn) if kwargs.get('writeCtf', True) and img.hasCTF(): ctfModelToRow(img.getCTF(), imgRow) # alignment is mandatory at this point, it shoud be check # and detected defaults if not passed at readSetOf.. level alignType = kwargs.get('alignType') if alignType != ALIGN_NONE: alignmentToRow(img.getTransform(), imgRow, alignType) if kwargs.get('writeAcquisition', True) and img.hasAcquisition(): acquisitionToRow(img.getAcquisition(), imgRow) # Write all extra labels to the row objectToRow(img, imgRow, {}, extraLabels=IMAGE_EXTRA_LABELS) # Provide a hook to be used if something is needed to be # done for special cases before converting image to row postprocessImageRow = kwargs.get('postprocessImageRow', None) if postprocessImageRow: postprocessImageRow(img, imgRow)
[docs]def rowToImage(imgRow, imgLabel, imgClass, **kwargs): """ Create an Image from a row of a metadata. """ img = imgClass() # Provide a hook to be used if something is needed to be # done for special cases before converting image to row preprocessImageRow = kwargs.get('preprocessImageRow', None) if preprocessImageRow: preprocessImageRow(img, imgRow) # Decompose Xmipp filename index, filename = xmippToLocation(imgRow.getValue(imgLabel)) img.setLocation(index, filename) if imgRow.containsLabel(emlib.MDL_REF): img.setClassId(imgRow.getValue(emlib.MDL_REF)) elif imgRow.containsLabel(emlib.MDL_REF3D): img.setClassId(imgRow.getValue(emlib.MDL_REF3D)) if kwargs.get('readCtf', True): img.setCTF(rowToCtfModel(imgRow)) # alignment is mandatory at this point, it shoud be check # and detected defaults if not passed at readSetOf.. level alignType = kwargs.get('alignType') if alignType != ALIGN_NONE: img.setTransform(rowToAlignment(imgRow, alignType)) if kwargs.get('readAcquisition', True): img.setAcquisition(rowToAcquisition(imgRow)) if kwargs.get('magnification', None): img.getAcquisition().setMagnification(kwargs.get("magnification")) setObjId(img, imgRow) # Read some extra labels rowToObject(imgRow, img, {}, extraLabels=IMAGE_EXTRA_LABELS + kwargs.get('extraLabels', [])) # Provide a hook to be used if something is needed to be # done for special cases before converting image to row postprocessImageRow = kwargs.get('postprocessImageRow', None) if postprocessImageRow: postprocessImageRow(img, imgRow) return img
[docs]def micrographToRow(mic, micRow, **kwargs): """ Set labels values from Micrograph mic to md row. """ imageToRow(mic, micRow, imgLabel=emlib.MDL_MICROGRAPH, **kwargs)
[docs]def rowToMicrograph(micRow, **kwargs): """ Create a Micrograph object from a row of Xmipp metadata. """ return rowToImage(micRow, emlib.MDL_MICROGRAPH, Micrograph, **kwargs)
[docs]def volumeToRow(vol, volRow, **kwargs): """ Set labels values from Micrograph mic to md row. """ imageToRow(vol, volRow, imgLabel=emlib.MDL_IMAGE, writeAcquisition=False, **kwargs)
[docs]def rowToVolume(volRow, **kwargs): """ Create a Volume object from a row of Xmipp metadata. """ return rowToImage(volRow, emlib.MDL_IMAGE, Volume, **kwargs)
[docs]def coordinateToRow(coord, coordRow, copyId=True): """ Set labels values from Coordinate coord to md row. """ if copyId: setRowId(coordRow, coord) objectToRow(coord, coordRow, COOR_DICT, extraLabels=COOR_EXTRA_LABELS) if coord.getMicId(): coordRow.setValue(emlib.MDL_MICROGRAPH, str(coord.getMicId()))
[docs]def rowToCoordinate(coordRow): """ Create a Coordinate from a row of a metadata. """ # Check that all required labels are present in the row if _containsAll(coordRow, COOR_DICT): coord = Coordinate() rowToObject(coordRow, coord, COOR_DICT, extraLabels=COOR_EXTRA_LABELS) # Setup the micId if is integer value try: coord.setMicId(int(coordRow.getValue(emlib.MDL_MICROGRAPH_ID))) except Exception: pass else: coord = None return coord
def _rowToParticle(partRow, particleClass, **kwargs): """ Create a Particle from a row of a metadata. """ # Since postprocessImage is intended to be after the object is # setup, we need to intercept it here and call it at the end postprocessImageRow = kwargs.get('postprocessImageRow', None) if postprocessImageRow: del kwargs['postprocessImageRow'] img = rowToImage(partRow, emlib.MDL_IMAGE, particleClass, **kwargs) img.setCoordinate(rowToCoordinate(partRow)) # copy micId if available # if not copy micrograph name if available try: if partRow.hasLabel(emlib.MDL_MICROGRAPH_ID): img.setMicId(partRow.getValue(emlib.MDL_MICROGRAPH_ID)) # elif partRow.hasLabel(emlib.MDL_MICROGRAPH): # micName = partRow.getValue(emlib.MDL_MICROGRAPH) # img._micrograph = micName # print("setting micname as %s" % micName) # img.printAll() # print("getAttributes1 %s" % img._micrograph) # print("getAttributes2 %s" % getattr(img, "_micrograph", 'kk') # else: # print("WARNING: No micname") except Exception as e: print("Warning:", e.message) if postprocessImageRow: postprocessImageRow(img, partRow) return img
[docs]def rowToParticle(partRow, **kwargs): return _rowToParticle(partRow, Particle, **kwargs)
[docs]def rowToMovieParticle(partRow, **kwargs): return _rowToParticle(partRow, MovieParticle, **kwargs)
[docs]def particleToRow(part, partRow, **kwargs): """ Set labels values from Particle to md row. """ imageToRow(part, partRow, emlib.MDL_IMAGE, **kwargs) coord = part.getCoordinate() if coord is not None: coordinateToRow(coord, partRow, copyId=False) if part.hasMicId(): partRow.setValue(emlib.MDL_MICROGRAPH_ID, int(part.getMicId())) partRow.setValue(emlib.MDL_MICROGRAPH, str(part.getMicId()))
[docs]def rowToClass(classRow, classItem): """ Method base to create a class2D, class3D or classVol from a row of a metadata """ setObjId(classItem, classRow, label=emlib.MDL_REF) if classRow.containsLabel(emlib.MDL_IMAGE): index, filename = xmippToLocation(classRow.getValue(emlib.MDL_IMAGE)) img = classItem.REP_TYPE() # class id should be set previously from MDL_REF # classItem.setObjId(classRow.getObjId()) # img.copyObjId(classItem) img.setLocation(index, filename) img.setSamplingRate(classItem.getSamplingRate()) classItem.setRepresentative(img) return classItem
[docs]def class2DToRow(class2D, classRow): """ Set labels values from Class2D to md row. """ if class2D.hasRepresentative(): index, filename = class2D.getRepresentative().getLocation() fn = locationToXmipp(index, filename) classRow.setValue(emlib.MDL_IMAGE, fn) n = int(len(class2D)) classRow.setValue(emlib.MDL_CLASS_COUNT, n) classRow.setValue(emlib.MDL_REF, int(class2D.getObjId())) classRow.setValue(emlib.MDL_ITEM_ID, int(class2D.getObjId()))
[docs]def ctfModelToRow(ctfModel, ctfRow): """ Set labels values from ctfModel to md row. """ # TODO: compatibility check remove eventually if ctfModel.hasResolution(): objectToRow(ctfModel, ctfRow, CTF_DICT, extraLabels=CTF_EXTRA_LABELS) if ctfModel.hasPhaseShift(): ctfRow.setValue(emlib.MDL_CTF_PHASE_SHIFT, ctfModel.getPhaseShift()) else: objectToRow(ctfModel, ctfRow, CTF_DICT_NORESOLUTION, extraLabels=CTF_EXTRA_LABELS)
[docs]def defocusGroupSetToRow(defocusGroup, defocusGroupRow): """ Set labels values from ctfModel to md row. """ objectToRow(defocusGroup, defocusGroupRow, CTF_DICT)
[docs]def setPsdFiles(ctfModel, ctfRow): """ Set the PSD files of CTF estimation related to this ctfModel. The values will be read from the ctfRow if present. """ for attr, label in CTF_PSD_DICT.items(): if ctfRow.containsLabel(label): setattr(ctfModel, attr, String(ctfRow.getValue(label)))
[docs]def rowToCtfModel(ctfRow): """ Create a CTFModel from a row of a metadata. """ # Check if the row has CTF values, this could be called from a xmipp # particles metadata if _containsAll(ctfRow, CTF_DICT_NORESOLUTION): # for compatibility reason ignore resolution and fitQuality # Instantiate Scipion CTF Model ctfModel = CTFModel() # Case for metadata coming with Xmipp resolution label # Populate Scipion CTF from metadata row (using mapping dictionary # plus extra labels if ctfRow.hasLabel(md.MDL_CTF_PHASE_SHIFT): ctfModel.setPhaseShift(ctfRow.getValue(md.MDL_CTF_PHASE_SHIFT, 0)) if ctfRow.containsLabel(emlib.label2Str(emlib.MDL_CTF_CRIT_MAXFREQ)): rowToObject(ctfRow, ctfModel, CTF_DICT, extraLabels=CTF_EXTRA_LABELS) else: rowToObject(ctfRow, ctfModel, CTF_DICT_NORESOLUTION, extraLabels=CTF_EXTRA_LABELS_PLUS_RESOLUTION) # Standarize defocus values ctfModel.standardize() # Set psd file names setPsdFiles(ctfModel, ctfRow) # ctfModel.setPhaseShift(0.0) # for consistency with ctfModel else: ctfModel = None return ctfModel
[docs]def acquisitionToRow(acquisition, ctfRow): """ Set labels values from acquisition to md row. """ objectToRow(acquisition, ctfRow, ACQUISITION_DICT)
[docs]def rowToAcquisition(acquisitionRow): """ Create an acquisition from a row of a metadata. """ if _containsAll(acquisitionRow, ACQUISITION_DICT): acquisition = Acquisition() rowToObject(acquisitionRow, acquisition, ACQUISITION_DICT) else: acquisition = None return acquisition
[docs]def readSetOfMicrographs(filename, micSet, **kwargs): readSetOfImages(filename, micSet, rowToMicrograph, **kwargs)
[docs]def writeSetOfMicrographs(micSet, filename, blockName='Micrographs', **kwargs): writeSetOfImages(micSet, filename, micrographToRow, blockName, **kwargs)
[docs]def readSetOfVolumes(filename, volSet, **kwargs): readSetOfImages(filename, volSet, rowToVolume, **kwargs)
[docs]def writeSetOfVolumes(volSet, filename, blockName='Volumes', **kwargs): writeSetOfImages(volSet, filename, volumeToRow, blockName, **kwargs)
[docs]def mdToCTFModel(mdIn, mic): ctfRow = rowFromMd(mdIn, mdIn.firstObject()) ctfObj = rowToCtfModel(ctfRow) setXmippAttributes(ctfObj, ctfRow, emlib.MDL_CTF_CRIT_NONASTIGMATICVALIDITY, emlib.MDL_CTF_CRIT_FIRSTMINIMUM_FIRSTZERO_DIFF_RATIO) ctfObj.setMicrograph(mic) return ctfObj
[docs]def readCTFModel(filename, mic): """ Read from Xmipp .ctfparam and create a CTFModel object. """ mdCtf = emlib.MetaData(filename) return mdToCTFModel(mdCtf, mic)
[docs]def openMd(fn, state='Manual'): # We are going to write metadata directy to file to do it faster f = open(fn, 'w') ismanual = state == 'Manual' block = 'data_particles' if ismanual else 'data_particles_auto' s = """# XMIPP_STAR_1 * # data_header loop_ _pickingMicrographState %s %s loop_ _itemId _enabled _xcoor _ycoor _cost _micrographId """ % (state, block) f.write(s) return f
[docs]def writeSetOfCoordinates(posDir, coordSet, ismanual=True, scale=1): state = 'Manual' if ismanual else 'Supervised' writeSetOfCoordinatesWithState(posDir, coordSet, state, scale)
[docs]def writeSetOfCoordinatesWithState(posDir, coordSet, state, scale=1): """ Write a pos file on metadata format for each micrograph on the coordSet. Params: posDir: the directory where the .pos files will be written. coordSet: the SetOfCoordinates that will be read. """ boxSize = coordSet.getBoxSize() or 100 # Create a dictionary with the pos filenames for each micrograph posDict = {} for mic in coordSet.iterMicrographs(): micIndex, micFileName = mic.getLocation() micName = os.path.basename(micFileName) if micIndex != NO_INDEX: micName = '%06d_at_%s' % (micIndex, micName) posFn = join(posDir, replaceBaseExt(micName, "pos")) posDict[mic.getObjId()] = posFn f = None lastMicId = None c = 0 for coord in coordSet.iterItems(orderBy='_micId'): micId = coord.getMicId() if micId != lastMicId: # we need to close previous opened file if f: f.close() c = 0 f = openMd(posDict[micId], state) lastMicId = micId c += 1 if scale != 1: x = coord.getX() * scale y = coord.getY() * scale else: x = coord.getX() y = coord.getY() f.write(" %06d 1 %d %d %d %06d\n" % (coord.getObjId(), x, y, 1, micId)) if f: f.close() # Write config.xmd metadata configFn = join(posDir, 'config.xmd') writeCoordsConfig(configFn, int(boxSize), state) return posDict.values()
[docs]def writeCoordsConfig(configFn, boxSize, state): """ Write the config.xmd file needed for Xmipp picker. Params: configFn: The filename were to store the configuration. boxSize: the box size in pixels for extraction. state: picker state """ # Write config.xmd metadata print("writeCoordsConfig: state=", state) mdCoord = emlib.MetaData() # Write properties block objId = mdCoord.addObject() mdCoord.setValue(emlib.MDL_PICKING_PARTICLE_SIZE, int(boxSize), objId) mdCoord.setValue(emlib.MDL_PICKING_STATE, state, objId) mdCoord.write('properties@%s' % configFn)
[docs]def writeMicCoordinates(mic, coordList, outputFn, isManual=True, getPosFunc=None): """ Write the pos file as expected by Xmipp with the coordinates of a given micrograph. Params: mic: input micrograph. coordList: list of (x, y) pairs of the mic coordinates. outputFn: output filename for the pos file . isManual: if the coordinates are 'Manual' or 'Supervised' getPosFunc: a function to get the positions from the coordinate, it can be useful for scaling the coordinates if needed. """ if getPosFunc is None: getPosFunc = lambda coord: coord.getPosition() state = 'Manual' if isManual else 'Supervised' f = openMd(outputFn, state) for coord in coordList: x, y = getPosFunc(coord) f.write(" %06d 1 %d %d %d %06d \n" % (coord.getObjId(), x, y, 1, mic.getObjId())) f.close()
[docs]def readSetOfCoordinates(outputDir, micSet, coordSet, readDiscarded=False, scale=1): """ Read from Xmipp .pos files. Params: outputDir: the directory where the .pos files are. It is also expected a file named: config.xmd in this directory where the box size can be read. micSet: the SetOfMicrographs to associate the .pos, which name should be the same of the micrographs. coordSet: the SetOfCoordinates that will be populated. readDiscarded: read only the coordinates with the MDL_ENABLE set at -1 scale: Factor to scale ONLY x,y coordinates, you are supposed to use an appropiate boxsize (you created the set) """ # Read the boxSize from the config.xmd metadata configfile = join(outputDir, 'config.xmd') if exists(configfile): mdCoords = emlib.MetaData('properties@' + join(outputDir, 'config.xmd')) boxSize = mdCoords.getValue(emlib.MDL_PICKING_PARTICLE_SIZE, mdCoords.firstObject()) coordSet.setBoxSize(int(boxSize)) #Only coordinates x,y are scaled, you are supposed to use an appropiate boxsize for mic in micSet: posFile = join(outputDir, replaceBaseExt(mic.getFileName(), 'pos')) readCoordinates(mic, posFile, coordSet, outputDir, readDiscarded, scale=scale) coordSet._xmippMd = String(outputDir)
[docs]def readCoordinates(mic, fileName, coordsSet, outputDir, readDiscarded=False, scale=1): posMd = readPosCoordinates(fileName, readDiscarded) # TODO: CHECK IF THIS LABEL IS STILL NECESSARY posMd.addLabel(md.MDL_ITEM_ID) for objId in posMd: # When do an union of two metadatas of coordinates and one of # them doesn't has MDL_ENABLED, the default vaule to is 0, # and its not allowed value. Maybe we need to solve this in xmipp # code. if posMd.getValue(md.MDL_ENABLED, objId) == 0: posMd.setValue(md.MDL_ENABLED, 1, objId) coord = rowToCoordinate(rowFromMd(posMd, objId)) coord.setMicrograph(mic) coord.setX(int(coord.getX()*scale)) coord.setY(int(coord.getY()*scale)) coordsSet.append(coord) posMd.setValue(md.MDL_ITEM_ID, int(coord.getObjId()), objId)
[docs]def readPosCoordinates(posFile, readDiscarded=False): """ Read the coordinates in .pos file and return corresponding metadata. There are two possible blocks with particles: particles: with manual/supervised particles particles_auto: with automatically picked particles. If posFile doesn't exist, the metadata will be empty readDiscarded: read only the coordinates in the particles_auto DB with the MDL_ENABLE set at -1 and a positive cost """ mData = md.MetaData() if exists(posFile): blocks = md.getBlocksInMetaDataFile(posFile) blocksToRead = ['particles_auto'] if readDiscarded \ else ['particles','particles_auto'] for b in blocksToRead: if b in blocks: mdAux = md.MetaData('%(b)s@%(posFile)s' % locals()) mData.unionAll(mdAux) if readDiscarded: for objId in mData: if mData.getValue(md.MDL_COST, objId)>0: enabled=mData.getValue(md.MDL_ENABLED, objId) mData.setValue(md.MDL_ENABLED, -1*enabled, objId) mData.removeDisabled() return mData
[docs]def readSetOfImages(filename, imgSet, rowToFunc, **kwargs): """read from Xmipp image metadata. filename: The metadata filename where the image are. imgSet: the SetOfParticles that will be populated. rowToFunc: this function will be used to convert the row to Object """ imgMd = emlib.MetaData(filename) # By default remove disabled items from metadata # be careful if you need to preserve the original number of items if kwargs.get('removeDisabled', True): imgMd.removeDisabled() # If the type of alignment is not sent through the kwargs # try to deduced from the metadata labels if 'alignType' not in kwargs: imgRow = rowFromMd(imgMd, imgMd.firstObject()) if _containsAny(imgRow, ALIGNMENT_DICT): if imgRow.containsLabel(emlib.MDL_ANGLE_TILT): kwargs['alignType'] = ALIGN_PROJ else: kwargs['alignType'] = ALIGN_2D else: kwargs['alignType'] = ALIGN_NONE if imgMd.size() > 0: for objId in imgMd: imgRow = rowFromMd(imgMd, objId) img = rowToFunc(imgRow, **kwargs) imgSet.append(img) imgSet.setHasCTF(img.hasCTF()) imgSet.setAlignment(kwargs['alignType'])
[docs]def setOfImagesToMd(imgSet, mdIn, imgToFunc, **kwargs): """ This function will fill Xmipp metadata from a SetOfMicrographs Params: imgSet: the set of images to be converted to metadata mdIn: metadata to be filled rowFunc: this function can be used to setup the row before adding to metadata. """ if 'alignType' not in kwargs: kwargs['alignType'] = imgSet.getAlignment() if 'where' in kwargs: where = kwargs['where'] for img in imgSet.iterItems(where=where): objId = mdIn.addObject() imgRow = md.Row() imgToFunc(img, imgRow, **kwargs) imgRow.writeToMd(mdIn, objId) else: for img in imgSet: objId = mdIn.addObject() imgRow = md.Row() imgToFunc(img, imgRow, **kwargs) imgRow.writeToMd(mdIn, objId)
[docs]def readAnglesFromMicrographs(micFile, anglesSet): """ Read the angles from a micrographs Metadata. """ micMd = emlib.MetaData(micFile) # micMd.removeDisabled() for objId in micMd: # TODO I have added the import from pyworkflow.em import Angles # because next line gave an error. I do not know if this # is the right Angle class angles = Angles() row = rowFromMd(micMd, objId) rowToObject(row, angles, ANGLES_DICT) angles.setObjId(micMd.getValue(emlib.MDL_ITEM_ID, objId)) anglesSet.append(angles)
[docs]def writeSetOfImages(imgSet, filename, imgToFunc, blockName='Images', **kwargs): """ This function will write a SetOfImages as a Xmipp metadata. Params: imgSet: the set of images to be written (particles, micrographs or volumes) filename: the filename where to write the metadata. rowFunc: this function can be used to setup the row before adding to metadata. """ mdSet = emlib.MetaData() setOfImagesToMd(imgSet, mdSet, imgToFunc, **kwargs) mdSet.write('%s@%s' % (blockName, filename))
[docs]def readSetOfParticles(filename, partSet, **kwargs): readSetOfImages(filename, partSet, rowToParticle, **kwargs)
[docs]def readSetOfMovieParticles(filename, partSet, **kwargs): readSetOfImages(filename, partSet, rowToMovieParticle, **kwargs)
[docs]def setOfParticlesToMd(imgSet, mdParts, **kwargs): setOfImagesToMd(imgSet, mdParts, particleToRow, **kwargs)
[docs]def setOfMicrographsToMd(imgSet, mdMics, **kwargs): setOfImagesToMd(imgSet, mdMics, micrographToRow, **kwargs)
[docs]def writeSetOfParticles(imgSet, filename, blockName='Particles', **kwargs): writeSetOfImages(imgSet, filename, particleToRow, blockName, **kwargs)
[docs]def writeCTFModel(ctfModel, ctfFile): """ Given a CTFModel object write as Xmipp ctfparam """ mdCtf = emlib.MetaData() objId = mdCtf.addObject() ctfRow = md.Row() ctfModelToRow(ctfModel, ctfRow) ctfRow.writeToMd(mdCtf, objId) mdCtf.setColumnFormat(False) mdCtf.write(ctfFile)
[docs]def writeSetOfCTFs(ctfSet, mdCTF): """ Write a ctfSet on metadata format. Params: ctfSet: the SetOfCTF that will be read. mdCTF: The file where metadata should be written. """ mdCtf = emlib.MetaData() for ctfModel in ctfSet: objId = mdCtf.addObject() ctfRow = md.Row() ctfRow.setValue(emlib.MDL_MICROGRAPH, ctfModel.getMicFile()) if ctfModel.getPsdFile(): ctfRow.setValue(emlib.MDL_PSD, ctfModel.getPsdFile()) ctfModelToRow(ctfModel, ctfRow) ctfRow.writeToMd(mdCtf, objId) mdCtf.write(mdCTF) ctfSet._xmippMd = String(mdCTF)
[docs]def writeSetOfDefocusGroups(defocusGroupSet, fnDefocusGroup): # also metadata """ Write a defocuGroupSet on metadata format. Params: defocusGroupSet: the SetOfDefocus that will be read. fnDefocusGroup: The file where defocusGroup should be written. """ mdDefGr = emlib.MetaData() for defocusGroup in defocusGroupSet: objId = mdDefGr.addObject() defocusGroupRow = md.Row() defocusGroupSetToRow(defocusGroup, defocusGroupRow) defocusGroupRow.setValue(emlib.MDL_CTF_GROUP, defocusGroup.getObjId()) defocusGroupRow.setValue(emlib.MDL_MIN, defocusGroup.getDefocusMin()) defocusGroupRow.setValue(emlib.MDL_MAX, defocusGroup.getDefocusMax()) defocusGroupRow.setValue(emlib.MDL_AVG, defocusGroup.getDefocusAvg()) defocusGroupRow.writeToMd(mdDefGr, objId) mdDefGr.write(fnDefocusGroup) defocusGroupSet._xmippMd = String(fnDefocusGroup)
[docs]def writeSetOfClasses2D(classes2DSet, filename, classesBlock='classes', writeParticles=True): """ This function will write a SetOfClasses2D as Xmipp metadata. Params: classes2DSet: the SetOfClasses2D instance. filename: the filename where to write the metadata. """ classFn = '%s@%s' % (classesBlock, filename) classMd = emlib.MetaData() classMd.write(classFn) # Empty write to ensure the classes is the first # block classRow = md.Row() for class2D in classes2DSet: class2DToRow(class2D, classRow) classRow.writeToMd(classMd, classMd.addObject()) if writeParticles: ref = class2D.getObjId() imagesFn = 'class%06d_images@%s' % (ref, filename) imagesMd = emlib.MetaData() imgRow = md.Row() if class2D.getSize() > 0: for img in class2D: particleToRow(img, imgRow) imgRow.writeToMd(imagesMd, imagesMd.addObject()) imagesMd.write(imagesFn, emlib.MD_APPEND) # Empty write to ensure the classes is the first block classMd.write(classFn, emlib.MD_APPEND)
[docs]def writeSetOfMicrographsPairs(uSet, tSet, filename): """ This function will write a MicrographsTiltPair as Xmipp metadata. Params: uSet: the untilted set of micrographs to be written tSet: the tilted set of micrographs to be written filename: the filename where to write the metadata. """ mdMics = emlib.MetaData() for micU, micT in izip(uSet, tSet): objId = mdMics.addObject() pairRow = md.Row() pairRow.setValue(emlib.MDL_ITEM_ID, int(micU.getObjId())) pairRow.setValue(emlib.MDL_MICROGRAPH, micU.getFileName()) pairRow.setValue(emlib.MDL_MICROGRAPH_TILTED, micT.getFileName()) pairRow.writeToMd(mdMics, objId) mdMics.write(filename)
def __readSetOfClasses(classBaseSet, readSetFunc, classesSet, filename, classesBlock='classes', **kwargs): """ Read a set of classes from an Xmipp metadata with the given convention of a block for classes and another block for each set of images assigned to each class. Params: classesSet: the SetOfClasses object that will be populated. filename: the file path where to read the Xmipp metadata. classesBlock (by default 'classes'): the block name of the classes group in the metadata. """ blocks = emlib.getBlocksInMetaDataFile(filename) classesMd = emlib.MetaData('%s@%s' % (classesBlock, filename)) # Provide a hook to be used if something is needed to be # done for special cases before converting row to class preprocessClass = kwargs.get('preprocessClass', None) postprocessClass = kwargs.get('postprocessClass', None) for objId in classesMd: classItem = classesSet.ITEM_TYPE() classRow = rowFromMd(classesMd, objId) # class id should be set in rowToClass function using MDL_REF classItem = rowToClass(classRow, classItem) classBaseSet.copyInfo(classItem, classesSet.getImages()) if preprocessClass: preprocessClass(classItem, classRow) classesSet.append(classItem) ref = classItem.getObjId() b = 'class%06d_images' % ref if b in blocks: readSetFunc('%s@%s' % (b, filename), classItem, **kwargs) if postprocessClass: postprocessClass(classItem, classRow) # Update with new properties of classItem such as _size classesSet.update(classItem)
[docs]def readSetOfClasses(classesSet, filename, classesBlock='classes', **kwargs): __readSetOfClasses(SetOfParticles, readSetOfParticles, classesSet, filename, classesBlock, **kwargs)
[docs]def readSetOfClasses2D(classes2DSet, filename, classesBlock='classes', **kwargs): """ Just a wrapper to readSetOfClasses. """ readSetOfClasses(classes2DSet, filename, classesBlock, **kwargs)
[docs]def readSetOfClasses3D(classes3DSet, filename, classesBlock='classes', **kwargs): """ Just a wrapper to readSetOfClasses. """ readSetOfClasses(classes3DSet, filename, classesBlock, **kwargs)
[docs]def writeSetOfClassesVol(classesVolSet, filename, classesBlock='classes'): """ This function will write a SetOfClassesVol as Xmipp metadata. Params: classesVolSet: the SetOfClassesVol instance. filename: the filename where to write the metadata. """ classFn = '%s@%s' % (classesBlock, filename) classMd = emlib.MetaData() classMd.write(classFn) # Empty write to ensure the classes is the first # block # FIXME: review implementation of this function since there are syntax # errors classRow = md.Row() for classVol in classesVolSet: # FIXME Where does this method come from classVolToRow(classVol, classRow) classRow.writeToMd(classMd, classMd.addObject()) ref = Class3D.getObjId() imagesFn = 'class%06d_images@%s' % (ref, filename) imagesMd = emlib.MetaData() imgRow = md.Row() for vol in classVol: volumeToRow(vol, imgRow) imgRow.writeToMd(imagesMd, imagesMd.addObject()) imagesMd.write(imagesFn, emlib.MD_APPEND) # Empty write to ensure the classes is the first block classMd.write(classFn, emlib.MD_APPEND)
[docs]def readSetOfClassesVol(classesVolSet, filename, classesBlock='classes', **kwargs): """read from Xmipp image metadata. fnImages: The metadata filename where the particles properties are. imgSet: the SetOfParticles that will be populated. hasCtf: is True if the ctf information exists. """ __readSetOfClasses(SetOfVolumes, readSetOfVolumes, classesVolSet, filename, classesBlock, **kwargs)
[docs]def writeSetOfMovies(moviesSet, filename, moviesBlock='movies'): """ This function will write a SetOfMovies as Xmipp metadata. Params: moviesSet: the SetOfMovies instance. filename: the filename where to write the metadata. """ for movie in moviesSet: ref = movie.getObjId() micrographsFn = 'movie%06d_micrographs@%s' % (ref, filename) micrographsMd = emlib.MetaData() micRow = md.Row() for mic in movie: micrographToRow(mic, micRow) micRow.writeToMd(micrographsMd, micrographsMd.addObject()) micrographsMd.write(micrographsFn, emlib.MD_APPEND)
[docs]def geometryFromMatrix(matrix, inverseTransform): from pwem.convert.transformations import (translation_from_matrix, euler_from_matrix) if inverseTransform: matrix = np.linalg.inv(matrix) shifts = -translation_from_matrix(matrix) else: shifts = translation_from_matrix(matrix) angles = -np.rad2deg(euler_from_matrix(matrix, axes='szyz')) return shifts, angles
[docs]def matrixFromGeometry(shifts, angles, inverseTransform): """ Create the transformation matrix from a given 2D shifts in X and Y...and the 3 euler angles. """ from pwem.convert import euler_matrix radAngles = -np.deg2rad(angles) M = euler_matrix(radAngles[0], radAngles[1], radAngles[2], 'szyz') if inverseTransform: M[:3, 3] = -shifts[:3] M = np.linalg.inv(M) else: M[:3, 3] = shifts[:3] return M
[docs]def rowToAlignment(alignmentRow, alignType): """ is2D == True-> matrix is 2D (2D images alignment) otherwise matrix is 3D (3D volume alignment or projection) invTransform == True -> for xmipp implies projection """ is2D = alignType == ALIGN_2D inverseTransform = alignType == ALIGN_PROJ if _containsAny(alignmentRow, ALIGNMENT_DICT): alignment = Transform() angles = np.zeros(3) shifts = np.zeros(3) flip = alignmentRow.getValue(emlib.MDL_FLIP) shifts[0] = alignmentRow.getValue(emlib.MDL_SHIFT_X, 0.) shifts[1] = alignmentRow.getValue(emlib.MDL_SHIFT_Y, 0.) if not is2D: angles[0] = alignmentRow.getValue(emlib.MDL_ANGLE_ROT, 0.) angles[1] = alignmentRow.getValue(emlib.MDL_ANGLE_TILT, 0.) shifts[2] = alignmentRow.getValue(emlib.MDL_SHIFT_Z, 0.) angles[2] = alignmentRow.getValue(emlib.MDL_ANGLE_PSI, 0.) if flip: angles[1] = angles[1]+180 # tilt + 180 angles[2] = - angles[2] # - psi, COSS: this is mirroring X shifts[0] = -shifts[0] # -x else: psi = alignmentRow.getValue(emlib.MDL_ANGLE_PSI, 0.) rot = alignmentRow.getValue(emlib.MDL_ANGLE_ROT, 0.) if rot != 0. and psi != 0: print("HORROR rot and psi are different from zero in 2D case") angles[0] = \ alignmentRow.getValue(emlib.MDL_ANGLE_PSI, 0.)\ + alignmentRow.getValue(emlib.MDL_ANGLE_ROT, 0.) matrix = matrixFromGeometry(shifts, angles, inverseTransform) if flip: if alignType == ALIGN_2D: matrix[0, :2] *= -1. # invert only the first two columns # keep x matrix[2, 2] = -1. # set 3D rot elif alignType == ALIGN_3D: matrix[0, :3] *= -1. # now, invert first line excluding x matrix[3, 3] *= -1. elif alignType == ALIGN_PROJ: pass alignment.setMatrix(matrix) # FIXME: now are also storing the alignment parameters since # the conversions to the Transform matrix have not been extensively # tested. # After this, we should only keep the matrix # for paramName, label in ALIGNMENT_DICT.iter(): # if alignmentRow.hasLabel(label): # setattr(alignment, paramName, # alignmentRow.getValueAsObject(label)) else: alignment = None return alignment
[docs]def alignmentToRow(alignment, alignmentRow, alignType): """ is2D == True-> matrix is 2D (2D images alignment) otherwise matrix is 3D (3D volume alignment or projection) invTransform == True -> for xmipp implies projection -> for xmipp implies alignment """ if alignment is None: return is2D = alignType == ALIGN_2D inverseTransform = alignType == ALIGN_PROJ # only flip is meaninfull if 2D case # in that case the 2x2 determinant is negative flip = False matrix = alignment.getMatrix() if alignType == ALIGN_2D: # get 2x2 matrix and check if negative flip = bool(np.linalg.det(matrix[0:2, 0:2]) < 0) if flip: matrix[0, :2] *= -1. # invert only the first two columns keep x matrix[2, 2] = 1. # set 3D rot else: pass elif alignType == ALIGN_3D: flip = bool(np.linalg.det(matrix[0:3, 0:3]) < 0) if flip: matrix[0, :4] *= -1. # now, invert first line including x matrix[3, 3] = 1. # set 3D rot else: pass else: flip = bool(np.linalg.det(matrix[0:3, 0:3]) < 0) if flip: raise Exception("the det of the transformation matrix is " "negative. This is not a valid transformation " "matrix for Scipion.") shifts, angles = geometryFromMatrix(matrix, inverseTransform) alignmentRow.setValue(emlib.MDL_SHIFT_X, shifts[0]) alignmentRow.setValue(emlib.MDL_SHIFT_Y, shifts[1]) if is2D: angle = angles[0] + angles[2] alignmentRow.setValue(emlib.MDL_ANGLE_PSI, angle) else: # if alignType == ALIGN_3D: alignmentRow.setValue(emlib.MDL_SHIFT_Z, shifts[2]) alignmentRow.setValue(emlib.MDL_ANGLE_ROT, angles[0]) alignmentRow.setValue(emlib.MDL_ANGLE_TILT, angles[1]) alignmentRow.setValue(emlib.MDL_ANGLE_PSI, angles[2]) alignmentRow.setValue(emlib.MDL_FLIP, flip)
[docs]def fillClasses(clsSet, updateClassCallback=None): """ Give an empty SetOfClasses (either 2D or 3D). Iterate over the input images and append to the corresponding class. It is important that each image has the classId properly set. """ clsDict = {} # Dictionary to store the (classId, classSet) pairs inputImages = clsSet.getImages() for img in inputImages: ref = img.getClassId() if ref is None: raise Exception('Particle classId is None!!!') # Register a new class set if the ref was not found. if ref not in clsDict: classItem = clsSet.ITEM_TYPE(objId=ref) rep = clsSet.REP_TYPE() classItem.setRepresentative(rep) clsDict[ref] = classItem classItem.copyInfo(inputImages) classItem.setAcquisition(inputImages.getAcquisition()) if updateClassCallback is not None: updateClassCallback(classItem) clsSet.append(classItem) else: # Try to get the class set given its ref number classItem = clsDict[ref] classItem.append(img) for classItem in clsDict.values(): clsSet.update(classItem)
# FIXME: USE THIS FUNCTION AS THE WAY TO CREATE CLASSES, REMOVE THE NEXT ONE
[docs]def createClassesFromImages2(inputImages, inputMd, classesFn, ClassType, classLabel, alignType, getClassFn=None, preprocessImageRow=None): """ From an intermediate X.xmd file produced by xmipp, create the set of classes in which those images are classified. Params: inputImages: the SetOfImages that were initially classified by relion. inputMd: the filename metadata. classesFn: filename where to write the classes. ClassType: the output type of the classes set ( usually SetOfClass2D or 3D ) classLabel: label that is the class reference. classFnTemplate: the template to get the classes averages filenames iter: the iteration number, just used in Class template """ mdIter = iterMdRows(inputMd) clsDict = {} # Dictionary to store the (classId, classSet) pairs clsSet = ClassType(filename=classesFn) clsSet.setImages(inputImages) hasCtf = inputImages.hasCTF() sampling = inputImages.getSamplingRate() for img, row in izip(inputImages, mdIter): ref = row.getValue(emlib.MDL_REF) if ref is None: raise Exception('MDL_REF not found in metadata: %s' % inputMd) if ref not in clsDict: # Register a new class set if the ref was not # found. classItem = clsSet.ITEM_TYPE(objId=ref) if getClassFn is None: refFn = '' else: refFn = getClassFn(ref) refLocation = xmippToLocation(refFn) rep = clsSet.REP_TYPE() rep.setLocation(refLocation) rep.setSamplingRate(sampling) classItem.setRepresentative(rep) clsDict[ref] = classItem classItem.copyInfo(inputImages) classItem.setAlignment(alignType) classItem.setAcquisition(inputImages.getAcquisition()) clsSet.append(classItem) else: # Try to get the class set given its ref number classItem = clsDict[ref] img.setTransform(rowToAlignment(row, alignType)) classItem.append(img) for classItem in clsDict.values(): clsSet.update(classItem) clsSet.write() return clsSet
# FIXME: remove this function and use previous one
[docs]def createClassesFromImages(inputImages, inputMd, classesFn, ClassType, classLabel, classFnTemplate, alignType, iter, preprocessImageRow=None): """ From an intermediate X.xmd file produced by xmipp, create the set of classes in which those images are classified. Params: inputImages: the SetOfImages that were initially classified by relion. inputMd: the filename metadata. classesFn: filename where to write the classes. ClassType: the output type of the classes set (usually SetOfClass2D or 3D) classLabel: label that is the class reference. classFnTemplate: the template to get the classes averages filenames iter: the iteration number, just used in Class template """ # We assume here that the volumes (classes3d) are in the same folder than # imgsFn # rootDir here is defined to be used expanding locals() if "@" in inputMd: inputFn = inputMd.split('@')[1] tmpDir = os.path.dirname(inputFn) else: tmpDir = os.path.dirname(inputMd) def getClassFn(ref): args = {'rootDir': tmpDir, 'iter': iter, 'ref': ref} return classFnTemplate % args createClassesFromImages2(inputImages, inputMd, classesFn, ClassType, classLabel, alignType, getClassFn=getClassFn, preprocessImageRow=preprocessImageRow)
[docs]def createItemMatrix(item, row, align): item.setTransform(rowToAlignment(row, alignType=align))
[docs]def readShiftsMovieAlignment(xmdFn): shiftsMd = md.MetaData(xmdFn) shiftsMd.removeDisabled() return (shiftsMd.getColumnValues(md.MDL_SHIFT_X), shiftsMd.getColumnValues(md.MDL_SHIFT_Y))
[docs]def writeShiftsMovieAlignment(movie, xmdFn, s0, sN): movieAlignment = movie.getAlignment() shiftListX, shiftListY = movieAlignment.getShifts() # Generating metadata for global shifts a0, aN = movieAlignment.getRange() globalShiftsMD = emlib.MetaData() alFrame = a0 if s0 < a0: for i in range(s0, a0): objId = globalShiftsMD.addObject() imgFn = locationToXmipp(i, getMovieFileName(movie)) globalShiftsMD.setValue(emlib.MDL_IMAGE, imgFn, objId) globalShiftsMD.setValue(emlib.MDL_SHIFT_X, 0.0, objId) globalShiftsMD.setValue(emlib.MDL_SHIFT_Y, 0.0, objId) for shiftX, shiftY in izip(shiftListX, shiftListY): if alFrame >= s0 and alFrame <= sN: objId = globalShiftsMD.addObject() imgFn = locationToXmipp(alFrame, getMovieFileName(movie)) globalShiftsMD.setValue(emlib.MDL_IMAGE, imgFn, objId) globalShiftsMD.setValue(emlib.MDL_SHIFT_X, shiftX, objId) globalShiftsMD.setValue(emlib.MDL_SHIFT_Y, shiftY, objId) alFrame += 1 if sN > aN: for j in range(aN, sN): objId = globalShiftsMD.addObject() imgFn = locationToXmipp(j+1, getMovieFileName(movie)) globalShiftsMD.setValue(emlib.MDL_IMAGE, imgFn, objId) globalShiftsMD.setValue(emlib.MDL_SHIFT_X, 0.0, objId) globalShiftsMD.setValue(emlib.MDL_SHIFT_Y, 0.0, objId) globalShiftsMD.write(xmdFn)
[docs]def writeMovieMd(movie, outXmd, f1, fN, useAlignment=False): movieMd = md.MetaData() frame = movie.clone() # get some info about the movie # problem is, that it can come from a movie set, and some # values might refer to different movie, e.g. no of frames :( firstFrame, _, frameIndex = movie.getFramesRange() # info from the movie set _, _ , lastFrame = movie.getDim() # actual no. of frame in the current movie lastFrame += 1 # (convert no. of frames to index, one-initiated) if lastFrame == 0: # this condition is for old SetOfMovies, that has lastFrame = 0. frames = movie.getNumberOfFrames() if frames is not None: lastFrame = frames if f1 < firstFrame or fN > lastFrame: raise Exception("Frame range could not be greater" " than the movie one.") ih = ImageHandler() if useAlignment: alignment = movie.getAlignment() if alignment is None: raise Exception("Can not write alignment for movie. ") a0, aN = alignment.getRange() if a0 < firstFrame or aN > lastFrame: raise Exception("Trying to write frames which have not been aligned.") shiftListX, shiftListY = alignment.getShifts() row = md.Row() stackIndex = frameIndex + (f1 - firstFrame) for i in range(f1, fN+1): frame.setIndex(stackIndex) row.setValue(md.MDL_IMAGE, ih.locationToXmipp(frame)) if useAlignment: shiftIndex = i - firstFrame row.setValue(emlib.MDL_SHIFT_X, shiftListX[shiftIndex]) row.setValue(emlib.MDL_SHIFT_Y, shiftListY[shiftIndex]) row.addToMd(movieMd) stackIndex += 1 movieMd.write(outXmd)
[docs]def createParamPhantomFile(filename, dimX, partSetMd, phFlip=False, ctfCorr=False): f = open(filename, 'w') str = "# XMIPP_STAR_1 *\n#\ndata_block1\n_dimensions2D '%d %d'\n" % \ (dimX, dimX) str += "_projAngleFile %s\n" % partSetMd str += "_ctfPhaseFlipped %d\n_ctfCorrected %d\n" % (phFlip, ctfCorr) str += "_applyShift 1\n_noisePixelLevel '0 0'\n" f.write(str) f.close()