Repository: I2PC/scipion-em-xmipp (language: python)

Issue:
Update CHANGES.md



Discussion / hints:

Fix patch:
diff --git a/.github/workflows/build.yml b/.github/workflows/build.yml
index 991b4c184..c32ead33d 100644
--- a/.github/workflows/build.yml
+++ b/.github/workflows/build.yml
@@ -11,6 +11,10 @@ jobs:
   build:
     # Specify the machine to run the job on
     runs-on: ubuntu-22.04
+    
+    #Avoid send statistics
+    env:
+      SEND_INSTALLATION_STATISTICS: "False"
 
     # Define the steps to be taken in the job
     steps:    
diff --git a/CHANGES.md b/CHANGES.md
index 85cb5a49c..d548bd3c3 100644
--- a/CHANGES.md
+++ b/CHANGES.md
@@ -1,3 +1,41 @@
+## Release 3.24.06 - Oceanus
+   - New protocols
+      - movie_alignment_consensus: protocol for the comparison of two global alignments
+      - PCA classification: 2D Classification method that can works in streaming and in static
+
+   - Protocols updated
+      - deep_center: Deep center calculates an approximate center for the particles.
+      - validate_fscq: Added cif support
+      - reconstruct_fourier: Enable reconstructing with halves in the reconstruct fourier
+      - simulate_ctf: Simulate CTF can add noise before and after CTF.
+      - extract_particles: Local average added
+      - compare_reprojections: Downsampling option,autodown sampling, allows input 2D in several formats, allows input 3D in several formats, ranking option available, extraction option available
+      - cl2d: ensuring a save classes creating and adding possible output
+   - Protocols fixed
+      - resolution_deepres:now works with queue system
+      - ProcessVolumes: Fix header for resize voxels
+      - resolution_bfactor:Fix local resolution to local bfactor  
+      - resolution_bfactor: create output for res bfac
+      - extract_particles: Bug fix in the downsampling factor of extract particles
+      - core_analysis: define correctly the particles pointer to have indirect pointer
+      - tilt_analysis: estimate automatically the window size new option
+      - movie_doseanalysis: update correctly the outputSet
+      - extract_particles: Check if resizing is enabled
+      - particle_pick_remove_duplicates: Remove duplicates coordinates 
+   - Protocols deprecated (For more details visit [this](https://github.com/I2PC/xmipp/wiki/Deprecating-programs-and-protocols))
+      - deep_denoissing
+      - metaprotocol_create_subset
+      - metaprotocol_golden_highres
+      - solid_angle
+      - split_volume
+   - More scipio-em-xmipp
+      - Solve Sonar Cloud reported bugs
+      - Flexible import of pyQT on chimera script
+      - Removing tkAgg matplotlib requirement in viewer_structure_map
+      - DLTK installation bug local variable and Fix use gpu
+      - updated scikit-learn version
+
+
 ## Release 3.23.11 - Nereus
    - New protocols
       - Volume local adjustment
diff --git a/xmipp3/__init__.py b/xmipp3/__init__.py
index 6ddaa2dfb..a0f16ff6b 100644
--- a/xmipp3/__init__.py
+++ b/xmipp3/__init__.py
@@ -41,8 +41,8 @@
 _logo = "xmipp_logo" + ("" if type_of_version == 'release' else '_devel') + '.png'
 
 _references = ['delaRosaTrevin2013', 'Sorzano2013', 'Strelak2021']
-_currentBinVersion = '3.23.11.0'
-__version__ = _currentBinVersion[2:] + ".1"  # Set this to ".0" on each xmipp binary release, otherwise increase it --> ".1", ".2", ...
+_currentBinVersion = '3.24.06.0'
+__version__ = _currentBinVersion[2:] + ".0"  # Set this to ".0" on each xmipp binary release, otherwise increase it --> ".1", ".2", ...
 
 # Requirement version variables
 NVIDIA_DRIVERS_MINIMUM_VERSION = 450
@@ -230,7 +230,6 @@ def installDeepLearningToolkit(plugin, env):
             msg = ("Tensorflow installed without GPU. Just CPU computations "
                    "enabled (slow computations).")
             cudaMsgs.append(msg)
-            useGpu = False
 
     if nvidiaDriverVer is not None:
         preMsgs.append("CUDA support found. Driver version: %s" % nvidiaDriverVer)
diff --git a/xmipp3/convert/convert.py b/xmipp3/convert/convert.py
index be597af8c..54fa31a1c 100644
--- a/xmipp3/convert/convert.py
+++ b/xmipp3/convert/convert.py
@@ -830,30 +830,35 @@ def writeSetOfCoordinatesWithState(posDir, coordSet, state, scale=1):
 
     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()
+    try:
+        for coord in coordSet.iterItems(orderBy='_micId'):
+            micId = coord.getMicId()
+
+            if micId != lastMicId:
+                # we need to close previous opened file
+                if f:
+                    f.close()
+                f = openMd(posDict[micId], state)
+                lastMicId = micId
+
+            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))
+
+    except (FileNotFoundError, PermissionError) as error:
+        print(f"Error occurred: {error}")
+    except Exception as e:
+        print(f"Unexpected error occurred: {e}")
+    finally:
+        if f:
+            f.close()
 
     # Write config.xmd metadata
     configFn = join(posDir, 'config.xmd')
@@ -1386,7 +1391,7 @@ def rowToAlignment(alignmentRow, alignType):
         else:
             psi = alignmentRow.getValue(emlib.MDL_ANGLE_PSI, 0.)
             rot = alignmentRow.getValue(emlib.MDL_ANGLE_ROT, 0.)
-            if rot != 0. and psi != 0:
+            if not np.isclose(rot, 0., atol=1e-6 ) and not np.isclose(psi, 0., atol=1e-6):
                 print("HORROR rot and psi are different from zero in 2D case")
             angles[0] = \
                 alignmentRow.getValue(emlib.MDL_ANGLE_PSI, 0.)\
diff --git a/xmipp3/protocols/protocol_deep_denoising.py b/xmipp3/legacy/protocols/protocol_deep_denoising.py
similarity index 100%
rename from xmipp3/protocols/protocol_deep_denoising.py
rename to xmipp3/legacy/protocols/protocol_deep_denoising.py
diff --git a/xmipp3/protocols/protocol_metaprotocol_create_subset.py b/xmipp3/legacy/protocols/protocol_metaprotocol_create_subset.py
similarity index 100%
rename from xmipp3/protocols/protocol_metaprotocol_create_subset.py
rename to xmipp3/legacy/protocols/protocol_metaprotocol_create_subset.py
diff --git a/xmipp3/protocols/protocol_metaprotocol_golden_highres.py b/xmipp3/legacy/protocols/protocol_metaprotocol_golden_highres.py
similarity index 100%
rename from xmipp3/protocols/protocol_metaprotocol_golden_highres.py
rename to xmipp3/legacy/protocols/protocol_metaprotocol_golden_highres.py
diff --git a/xmipp3/protocols/protocol_solid_angles.py b/xmipp3/legacy/protocols/protocol_solid_angles.py
similarity index 100%
rename from xmipp3/protocols/protocol_solid_angles.py
rename to xmipp3/legacy/protocols/protocol_solid_angles.py
diff --git a/xmipp3/protocols/protocol_split_volume.py b/xmipp3/legacy/protocols/protocol_split_volume.py
similarity index 100%
rename from xmipp3/protocols/protocol_split_volume.py
rename to xmipp3/legacy/protocols/protocol_split_volume.py
diff --git a/xmipp3/viewers/viewer_metaprotocol_golden_highres.py b/xmipp3/legacy/viewers/viewer_metaprotocol_golden_highres.py
similarity index 100%
rename from xmipp3/viewers/viewer_metaprotocol_golden_highres.py
rename to xmipp3/legacy/viewers/viewer_metaprotocol_golden_highres.py
diff --git a/xmipp3/viewers/viewer_solid_angles.py b/xmipp3/legacy/viewers/viewer_solid_angles.py
similarity index 100%
rename from xmipp3/viewers/viewer_solid_angles.py
rename to xmipp3/legacy/viewers/viewer_solid_angles.py
diff --git a/xmipp3/viewers/viewer_split_volume.py b/xmipp3/legacy/viewers/viewer_split_volume.py
similarity index 100%
rename from xmipp3/viewers/viewer_split_volume.py
rename to xmipp3/legacy/viewers/viewer_split_volume.py
diff --git a/xmipp3/protocols.conf b/xmipp3/protocols.conf
index 1fd80b05b..f908d4a76 100644
--- a/xmipp3/protocols.conf
+++ b/xmipp3/protocols.conf
@@ -3,6 +3,7 @@ Protocols SPA = [
 	{"tag": "section", "text": "Movies", "openItem": "False", "children": [
 		{"tag": "protocol", "value": "XmippProtFlexAlign", "text": "default"},
 		{"tag": "protocol", "value": "XmippProtOFAlignment", "text": "default"},
+		{"tag": "protocol", "value": "XmippProtConsensusMovieAlignment", "text": "default"},
 		{"tag": "protocol", "value": "XmippProtMovieMaxShift", "text": "default"},
 		{"tag": "protocol", "value": "XmippProtSplitFrames", "text": "default"}	,
 		{"tag": "protocol", "value": "XmippProtMovieGain", "text": "default"},
@@ -56,13 +57,15 @@ Protocols SPA = [
 	]},
 	{"tag": "section", "text": "2D", "children": [
 		{"tag": "protocol_group", "text": "Align", "openItem": "False", "children": [
-			{"tag": "protocol", "value": "XmippProtCL2DAlign",   "text": "default"},		
+			{"tag": "protocol", "value": "XmippProtCL2DAlign",   "text": "default"},			
+			{"tag": "protocol", "value": "XmippProtDeepCenter",   "text": "default"},		
 			{"tag": "section", "text": "more", "openItem": "False", "children": [
 				{"tag": "protocol", "value": "ProtAlignmentAssign", "text": "default"}
 			]}
 		]},
 		{"tag": "protocol_group", "text": "Classify", "openItem": "False", "children": [
-			{"tag": "protocol", "value": "XmippProtCL2D",   "text": "default"},			
+			{"tag": "protocol", "value": "XmippProtCL2D",   "text": "default"},	
+			{"tag": "protocol", "value": "XmippProtClassifyPca",   "text": "default"},		
 			{"tag": "protocol", "value": "XmippProtCoreAnalysis",   "text": "default"},
 			{"tag": "section", "text": "more", "openItem": "False", "children": [
 				{"tag": "protocol", "value": "XmippProtML2D", "text": "default"},
diff --git a/xmipp3/protocols/__init__.py b/xmipp3/protocols/__init__.py
index 9c8244e0b..7cb2acd55 100644
--- a/xmipp3/protocols/__init__.py
+++ b/xmipp3/protocols/__init__.py
@@ -39,6 +39,8 @@
 from .protocol_break_symmetry import XmippProtAngBreakSymmetry
 from .protocol_cl2d_align import XmippProtCL2DAlign
 from .protocol_cl2d import XmippProtCL2D
+from .protocol_classify_pca import XmippProtClassifyPca
+from .protocol_classify_pca_streaming import XmippProtClassifyPcaStreaming
 #from .protocol_classify_kmeans2d import XmippProtKmeansClassif2D
 from .protocol_ctf_defocus_group import XmippProtCTFDefocusGroup
 from .protocol_compare_reprojections import XmippProtCompareReprojections
@@ -51,7 +53,7 @@
 from .protocol_ctf_correct_wiener2d import XmippProtCTFCorrectWiener2D
 from .protocol_consensus_classes import XmippProtConsensusClasses
 from .protocol_denoise_particles import XmippProtDenoiseParticles
-from .protocol_deep_micrograph_screen import  XmippProtDeepMicrographScreen
+from .protocol_deep_micrograph_screen import XmippProtDeepMicrographScreen
 from .protocol_eliminate_empty_images import (XmippProtEliminateEmptyParticles,
                                               XmippProtEliminateEmptyClasses)
 from .protocol_enrich import XmippProtEnrich
@@ -63,6 +65,7 @@
 from .protocol_kerdensom import XmippProtKerdensom
 from .protocol_ml2d import XmippProtML2D
 from .protocol_movie_gain import XmippProtMovieGain
+from .protocol_movie_alignment_consensus import XmippProtConsensusMovieAlignment
 from .protocol_flexalign import XmippProtFlexAlign
 from .protocol_movie_opticalflow import XmippProtOFAlignment, ProtMovieAlignment
 from .protocol_movie_max_shift import XmippProtMovieMaxShift
@@ -103,8 +106,6 @@
 from .protocol_shift_particles import XmippProtShiftParticles
 from .protocol_shift_volume import XmippProtShiftVolume
 from .protocol_simulate_ctf import XmippProtSimulateCTF
-from .protocol_solid_angles import XmippProtSolidAngles
-from .protocol_split_volume import XmippProtSplitvolume
 from .protocol_subtract_projection import XmippProtSubtractProjection
 from .protocol_subtract_projection import XmippProtBoostParticles
 from .protocol_tilt_analysis import XmippProtTiltAnalysis
@@ -117,8 +118,6 @@
 from .protocol_write_testC import XmippProtWriteTestC
 from .protocol_write_testP import XmippProtWriteTestP
 from .protocol_generate_reprojections import XmippProtGenerateReprojections
-from .protocol_deep_denoising import XmippProtDeepDenoising  # deepDenoising has to be after XmippProtGenerateReprojections and XmippProtCompareReprojections
-from .protocol_metaprotocol_create_subset import XmippMetaProtCreateSubset
 from .protocol_volume_deform_zernike3d import XmippProtVolumeDeformZernike3D
 from .protocol_structure_map_zernike3d import XmippProtStructureMapZernike3D
 from .protocol_align_volume_and_particles import XmippProtAlignVolumeParticles
@@ -126,7 +125,6 @@
 from .protocol_analyze_local_ctf import XmippProtAnalyzeLocalCTF
 from .protocol_consensus_local_ctf import XmippProtConsensusLocalCTF
 from .protocol_particle_pick_remove_duplicates import XmippProtPickingRemoveDuplicates
-from .protocol_metaprotocol_golden_highres import XmippMetaProtGoldenHighRes
 # from .protocol_apply_deformation_zernike3d import XmippProtApplyZernike3D
 # from .protocol_kmeans_clustering import XmippProtKmeansSPH
 from .protocol_structure_map import XmippProtStructureMap
@@ -136,4 +134,4 @@
 from .protocol_volume_local_adjust import XmippProtLocalVolAdj
 from .protocol_classes_2d_mapping import XmippProtCL2DMap
 from .protocol_deep_hand import XmippProtDeepHand
-from .protocol_deep_center_assignment import XmippProtDeepCenter, XmippProtDeepGlobalAssignment
+from .protocol_deep_center import XmippProtDeepCenter
diff --git a/xmipp3/protocols/protocol_angular_graph_consistency.py b/xmipp3/protocols/protocol_angular_graph_consistency.py
index 9d2020961..90b4cd6e1 100644
--- a/xmipp3/protocols/protocol_angular_graph_consistency.py
+++ b/xmipp3/protocols/protocol_angular_graph_consistency.py
@@ -521,7 +521,7 @@ def _summary(self):
         if (not hasattr(self, 'outputParticles')):
             summary.append("Output not ready yet.")
         else:
-            if self.percentage == 100.0:
+            if self.percentage == 100:
                 text = 'After validation, all of the images are likely to be within the realiable assignment zone'
             else:
                 text = 'After validation, a %.2f' % self.percentage
diff --git a/xmipp3/protocols/protocol_cl2d.py b/xmipp3/protocols/protocol_cl2d.py
index bf6c32c76..20c69da9f 100644
--- a/xmipp3/protocols/protocol_cl2d.py
+++ b/xmipp3/protocols/protocol_cl2d.py
@@ -36,6 +36,7 @@
 from pwem.protocols import ProtClassify2D
 from pwem.objects import SetOfClasses2D
 from pwem.constants import ALIGN_NONE, ALIGN_2D
+from pyworkflow import BETA, UPDATED, NEW, PROD
 
 
 from xmipp3.convert import (writeSetOfParticles, createItemMatrix,
@@ -52,6 +53,7 @@
 CL_ROBUST = 1
 
 # Classes keys
+OUTPUTCLASSES = 'outputClasses'
 CLASSES = ''
 CLASSES_CORE = '_core'
 CLASSES_STABLE_CORE = '_stable_core'
@@ -65,6 +67,11 @@ class XmippProtCL2D(ProtClassify2D):
     the original dataset into a given number of classes. """
     
     _label = 'cl2d'
+    _devStatus = UPDATED
+
+    _possibleOutputs = {OUTPUTCLASSES: SetOfClasses2D,
+                        OUTPUTCLASSES+CLASSES_CORE: SetOfClasses2D,
+                        OUTPUTCLASSES+CLASSES_STABLE_CORE: SetOfClasses2D}
     
     def __init__(self, **args):
         ProtClassify2D.__init__(self, **args)
@@ -195,7 +202,7 @@ def _insertAllSteps(self):
         else:
             initialClassesId = None
         
-        self._insertFunctionStep('convertInputStep',
+        self._insertFunctionStep(self.convertInputStep,
                                  self.inputParticles.get().getObjId(),
                                  initialClassesId)
 
@@ -229,7 +236,7 @@ def _insertAllSteps(self):
             args = self._defArgsCoreAnalisys()
             self._insertClassifySteps(program, args, subset=CLASSES_CORE)
             if self.analyzeRejected:
-                self._insertFunctionStep('analyzeOutOfCores', CLASSES_CORE)
+                self._insertFunctionStep(self.analyzeOutOfCores, CLASSES_CORE)
 
             if (self.numberOfClasses > (2 * self.numberOfInitialClasses.get())
                 and self.doStableCore): # Number of levels should be > 2
@@ -239,7 +246,7 @@ def _insertAllSteps(self):
                 self._insertClassifySteps(program, args,
                                           subset=CLASSES_STABLE_CORE)
                 if self.analyzeRejected:
-                    self._insertFunctionStep('analyzeOutOfCores',
+                    self._insertFunctionStep(self.analyzeOutOfCores,
                                              CLASSES_STABLE_CORE)
 
     def _insertClassifySteps(self, program, args, subset=CLASSES):
@@ -250,9 +257,9 @@ def _insertClassifySteps(self, program, args, subset=CLASSES):
         4. And create output
         """
         self._insertRunJobStep(program, args % self._params)
-        self._insertFunctionStep('evaluateClassesStep', subset)
-        self._insertFunctionStep('sortClassesStep', subset)
-        self._insertFunctionStep('createOutputStep', subset)
+        self._insertFunctionStep(self.evaluateClassesStep, subset)
+        self._insertFunctionStep(self.sortClassesStep, subset)
+        self._insertFunctionStep(self.createOutputStep, subset)
     
     #--------------------------- STEPS functions -------------------------------
     def convertInputStep(self, particlesId, classesId):
@@ -309,16 +316,15 @@ def createOutputStep(self, subset=''):
         resulting from the protocol execution.
         """
        
-        inputParticles = self.inputParticles.get()
         level = self._lastLevel()
         
         subsetFn = self._getFileName("level_classes", level=level, sub=subset)
         
         if exists(subsetFn):
-            classes2DSet = self._createSetOfClasses2D(inputParticles, subset)
+            classes2DSet = self._createSetOfClasses2D(self.inputParticles, subset)
             self._fillClassesFromLevel(classes2DSet, "last", subset)
     
-            result = {'outputClasses' + subset: classes2DSet}
+            result = {OUTPUTCLASSES + subset: classes2DSet}
             self._defineOutputs(**result)
             self._defineSourceRelation(self.inputParticles, classes2DSet)
     
diff --git a/xmipp3/protocols/protocol_classify_pca.py b/xmipp3/protocols/protocol_classify_pca.py
new file mode 100644
index 000000000..977d2f0d4
--- /dev/null
+++ b/xmipp3/protocols/protocol_classify_pca.py
@@ -0,0 +1,466 @@
+# ******************************************************************************
+# *
+# * Authors:     Erney Ramirez Aportela (eramirez@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'
+# *
+# ******************************************************************************
+import os
+import emtable
+import enum
+import numpy as np
+
+from pyworkflow import VERSION_3_0
+from pyworkflow.protocol.params import (PointerParam, StringParam, FloatParam,
+                                        BooleanParam, EnumParam, IntParam, GPU_LIST)
+from pyworkflow.protocol.constants import LEVEL_ADVANCED
+from pyworkflow.constants import BETA
+
+from pwem.protocols import ProtClassify2D
+from pwem.objects import SetOfClasses2D, Transform
+from pwem.constants import ALIGN_NONE, ALIGN_2D, ALIGN_PROJ, ALIGN_3D
+
+from xmipp3 import XmippProtocol
+from xmipp3.convert import (writeSetOfParticles, writeSetOfClasses2D, xmippToLocation,
+                            readSetOfParticles, matrixFromGeometry)
+
+
+class XMIPPCOLUMNS(enum.Enum):
+    # PARTICLES CONSTANTS
+    ctfVoltage = "ctfVoltage"  # 1
+    ctfDefocusU = "ctfDefocusU"  # 2
+    ctfDefocusV = "ctfDefocusV"  # 3
+    ctfDefocusAngle = "ctfDefocusAngle"  # 4
+    ctfSphericalAberration = "ctfSphericalAberration"  # 5
+    ctfQ0 = "ctfQ0"  # 6
+    ctfCritMaxFreq = "ctfCritMaxFreq"  # 7
+    ctfCritFitting = "ctfCritFitting"  # 8
+    enabled = "enabled"  # 9
+    image = "image"  # 10
+    itemId = "itemId"  # 11
+    micrograph = "micrograph"  # 12
+    micrographId = "micrographId"  # 13
+    scoreByVariance = "scoreByVariance"  # 14
+    scoreByGiniCoeff = "scoreByGiniCoeff"  # 15
+    xcoor = "xcoor"  # 16
+    ycoor = "ycoor"  # 17
+    ref = "ref"  # 18
+    anglePsi = "anglePsi"  # 19
+    angleRot = "angleRot"  # 20
+    angleTilt = "angleTilt"  # 21
+    shiftX = "shiftX"  # 22
+    shiftY = "shiftY"  # 23
+    shiftZ = "shiftZ"  # 24
+    flip = "flip"
+
+    # CLASSES CONSTANTS
+    classCount = "classCount"  # 3
+
+
+ALIGNMENT_DICT = {"shiftX": XMIPPCOLUMNS.shiftX.value,
+                  "shiftY": XMIPPCOLUMNS.shiftY.value,
+                  "shiftZ": XMIPPCOLUMNS.shiftZ.value,
+                  "flip": XMIPPCOLUMNS.flip.value,
+                  "anglePsi": XMIPPCOLUMNS.anglePsi.value,
+                  "angleRot": XMIPPCOLUMNS.angleRot.value,
+                  "angleTilt": XMIPPCOLUMNS.angleTilt.value
+                  }
+
+
+
+def updateEnviron(gpuNum):
+    """ Create the needed environment for pytorch programs. """
+    print("updating environ to select gpu %s" % (gpuNum))
+    if gpuNum == '':
+        os.environ['CUDA_VISIBLE_DEVICES'] = '0'
+    else:
+        os.environ['CUDA_VISIBLE_DEVICES'] = str(gpuNum)
+
+CONTRAST_AVERAGES_FILE = 'classes_contrast_classes.star'
+AVERAGES_IMAGES_FILE = 'classes_images.star'
+        
+class XmippProtClassifyPca(ProtClassify2D, XmippProtocol):
+    """ Classifies a set of images. """
+    
+    _label = '2D classification pca'
+    _lastUpdateVersion = VERSION_3_0
+    _conda_env = 'xmipp_pyTorch'
+    _devStatus = BETA
+    
+    # Mode
+    CREATE_CLASSES = 0
+    UPDATE_CLASSES = 1
+    
+    def __init__(self, **args):
+        ProtClassify2D.__init__(self, **args)
+
+    # --------------------------- DEFINE param functions ------------------------
+    def _defineParams(self, form):
+        form = self._defineCommonParams(form)
+        form.addParallelSection(threads=1, mpi=8)
+
+    def _defineCommonParams(self, form):
+        form.addHidden(GPU_LIST, StringParam, default='0',
+                       label="Choose GPU ID",
+                       help="GPU may have several cores. Set it to zero"
+                            " if you do not know what we are talking about."
+                            " First core index is 0, second 1 and so on.")
+
+        form.addSection(label='Input')
+
+        form.addParam('inputParticles', PointerParam,
+                      label="Input images",
+                      important=True, pointerClass='SetOfParticles',
+                      help='Select the input images to be classified.')
+        form.addParam('mode', EnumParam, choices=['create_classes', 'update_classes'],
+                      label="Create or update 2D classes?", default=self.CREATE_CLASSES,
+                      display=EnumParam.DISPLAY_HLIST,
+                      help='This option allows for the classification '
+                           'or simply alignment of particles into previously created classes.')
+        form.addParam('numberOfClasses', IntParam, default=50,
+                      condition="not mode",
+                      label='Number of classes:',
+                      help='Number of classes (or references) to be generated.')
+        form.addParam('initialClasses', PointerParam,
+                      label="Initial classes",
+                      condition="mode",
+                      pointerClass='SetOfClasses2D, SetOfAverages',
+                      help='Set of initial classes to start the classification')
+        form.addParam('correctCtf', BooleanParam, default=True, expertLevel=LEVEL_ADVANCED,
+                      label='Correct CTF?',
+                      help='If you set to *Yes*, the CTF of the experimental particles will be corrected')
+        form.addParam('mask', BooleanParam, default=True, expertLevel=LEVEL_ADVANCED,
+                      label='Use Gaussian Mask?',
+                      help='If you set to *Yes*, a gaussian mask is applied to the images.')
+        form.addParam('sigma', IntParam, default=-1, expertLevel=LEVEL_ADVANCED,
+                      label='sigma:', condition="mask",
+                      help='Sigma is the parameter that controls the dispersion or "width" of the curve.'
+                           ' If the parameter is set to -1, sigma = dim/3.')
+
+        form.addSection(label='Pca training')
+        form.addParam('resolution', FloatParam, label="max resolution", default=8,
+                      help='Maximum resolution to be consider for alignment')
+        form.addParam('coef', FloatParam, label="% variance", default=0.75, expertLevel=LEVEL_ADVANCED,
+                      help='Percentage of variance to determine the number of PCA components (between 0-1).'
+                           ' The higher the percentage, the higher the accuracy, but the calculation time increases.')
+        form.addParam('training', IntParam, default=40000,
+                      label="particles for training",
+                      help='Number of particles for PCA training')
+
+        return form
+
+    # --------------------------- INSERT steps functions ------------------------
+    def _insertAllSteps(self):
+
+        """ Mainly prepare the command line for call classification program"""
+    
+        updateEnviron(self.gpuList.get())
+        
+        self.imgsOrigXmd = self._getExtraPath('images_original.xmd')
+        self.imgsXmd = self._getTmpPath('images.xmd')
+        self.imgsFn = self._getTmpPath('images.mrc')
+        self.refXmd = self._getTmpPath('references.xmd')
+        self.ref = self._getTmpPath('references.mrcs')
+        self.sampling = self.inputParticles.get().getSamplingRate()
+        self.acquisition = self.inputParticles.get().getAcquisition()
+        mask = self.mask.get()
+        sigma = self.sigma.get()
+        if sigma == -1:
+            sigma = self.inputParticles.get().getDimensions()[0]/3
+        particlesTrain = min(self.training.get(), self.inputParticles.get().getSize())
+        resolution = self.resolution.get()
+        if resolution < 2 * self.sampling:
+            resolution = (2 * self.sampling) + 0.5
+
+    
+        self._insertFunctionStep('convertInputStep', 
+                                self.inputParticles.get(), self.imgsOrigXmd, self.imgsFn) #convert
+        
+        self._insertFunctionStep("pcaTraining", self.imgsFn, resolution, particlesTrain)
+        
+        self._insertFunctionStep("classification", self.imgsFn, self.numberOfClasses.get(), self.imgsOrigXmd, mask, sigma )
+    
+        self._insertFunctionStep('createOutputStep')
+        
+
+    #--------------------------- STEPS functions -------------------------------
+    def convertInputStep(self, input, outputOrig, outputMRC):
+        writeSetOfParticles(input, outputOrig)
+        
+        if self.mode == self.UPDATE_CLASSES: 
+            
+            if isinstance(self.initialClasses.get(), SetOfClasses2D):
+                writeSetOfClasses2D(self.initialClasses.get(),
+                                    self.refXmd, writeParticles=False)
+            else:
+                writeSetOfParticles(self.initialClasses.get(),
+                                    self.refXmd)
+                       
+            args = ' -i  %s -o %s  '%(self.refXmd, self.ref)
+            self.runJob("xmipp_image_convert", args, numberOfMpi=1)
+            
+        
+        if self.correctCtf: 
+            args = ' -i  %s -o %s --sampling_rate %s '%(outputOrig, outputMRC, self.sampling)
+            self.runJob("xmipp_ctf_correct_wiener2d", args, numberOfMpi=self.numberOfMpi.get())
+            
+        else:      
+            args = ' -i  %s -o %s  '%(outputOrig, outputMRC)
+            self.runJob("xmipp_image_convert", args, numberOfMpi=1) 
+        
+        
+    def pcaTraining(self, inputIm, resolutionTrain, numTrain):
+        args = ' -i %s  -s %s -hr %s -lr 530 -p %s -t %s -o %s/train_pca  --batchPCA'% \
+                (inputIm, self.sampling, resolutionTrain, self.coef.get(), numTrain, self._getExtraPath())
+
+        env = self.getCondaEnv()
+        env['LD_LIBRARY_PATH'] = ''
+        self.runJob("xmipp_classify_pca_train", args, numberOfMpi=1, env=env)
+        
+        
+    def classification(self, inputIm, numClass, stfile, mask, sigma):
+        args = ' -i %s -c %s -b %s/train_pca_bands.pt -v %s/train_pca_vecs.pt -o %s/classes -stExp %s' % \
+                (inputIm, numClass, self._getExtraPath(), self._getExtraPath(),  self._getExtraPath(),
+                 stfile)
+        if mask:
+            args += ' --mask --sigma %s '%(sigma) 
+            
+        if self.mode == self.UPDATE_CLASSES:
+            args += ' -r %s '%(self.ref)
+
+        env = self.getCondaEnv()
+        env['LD_LIBRARY_PATH'] = ''
+        self.runJob("xmipp_classify_pca", args, numberOfMpi=1, env=env)
+        
+        
+    def createOutputStep(self):
+        """ Store the SetOfClasses2D object
+        resulting from the protocol execution.
+        """
+    
+        inputParticles = self.inputParticles#.get()
+
+        classes2DSet = self._createSetOfClasses2D(inputParticles)
+        self._fillClassesFromLevel(classes2DSet)
+
+        self._defineOutputs(outputClasses=classes2DSet)
+        self._defineSourceRelation(self.inputParticles, classes2DSet)
+        
+        self.createOutputAverages(classes2DSet)
+        
+    def createOutputAverages(self, outputClasses):
+        classes = self._getExtraPath('classes.mrcs')
+        outRefs = self._createSetOfAverages()
+        outRefs.copyInfo(self.inputParticles.get())
+        outRefs.setSamplingRate(self.sampling)
+        readSetOfParticles(classes, outRefs)
+
+        outRefs.setAlignment(ALIGN_2D)
+        self._defineOutputs(outputAverages=outRefs)
+        self._defineSourceRelation(self.inputParticles, outRefs)    
+    
+    # --------------------------- INFO functions --------------------------------
+    def _validate(self):
+        """ Check if the installation of this protocol is correct.
+        Can't rely on package function since this is a "multi package" package
+        Returning an empty list means that the installation is correct
+        and there are not errors. If some errors are found, a list with
+        the error messages will be returned.
+        """
+        
+        errors = []
+        if self.inputParticles.get().getDimensions()[0] > 256:
+            errors.append("You should resize the particles."
+                          " Sizes smaller than 128 pixels are recommended.")
+        er = self.validateDLtoolkit()
+        if er:
+            errors.append(er)
+        return errors
+
+    def _warnings(self):
+        validateMsgs = []
+        if self.inputParticles.get().getDimensions()[0] > 128:
+            validateMsgs.append("Particle sizes equal to or less"
+                                " than 128 pixels are recommended.")
+        if self.inputParticles.get().getDimensions()[0] > 256:
+            validateMsgs.append("Particle sizes bigger than 256 may"
+                                " saturate the GPU memory.")
+        return validateMsgs
+
+    def _summary(self):
+        summary = []
+
+        if not hasattr(self, 'outputClasses'):
+            summary.append("Output classes not ready yet.")
+        else:
+            summary.append('Two output sets are obtained:')
+            summary.append('- The first one corresponds to the classes and should be used to select '
+                           ' the particles corresponding to each class. In this case, the classes are '
+                           ' displayed applying a contrast variation.')
+            summary.append('- The second one corresponds solely to the representative classes (outputAverages)')
+        return summary
+
+    # #--------------------------- UTILS functions -------------------------------
+    # EMTABLE IMPLEMENTATION
+    def _updateParticle(self, item, row):
+        if row is None:
+            self.info('Row is none finish updating particle')
+            setattr(item, "_appendItem", False)
+        else:
+            if item.getObjId() == row.get(XMIPPCOLUMNS.itemId.value):
+                item.setClassId(row.get(XMIPPCOLUMNS.ref.value))
+                item.setTransform(rowToAlignmentEmtable(row, ALIGN_2D))
+            else:
+                self.error('The particles ids are not synchronized')
+                setattr(item, "_appendItem", False)
+
+    def _updateClass(self, item):
+        classId = item.getObjId()
+        if classId in self._classesInfo:
+            index, fn, _ = self._classesInfo[classId]
+            item.setAlignment2D()
+            rep = item.getRepresentative()
+            rep.setLocation(index, fn)
+            rep.setSamplingRate(self.inputParticles.get().getSamplingRate())
+
+    def _createModelFile(self):
+        with open(self._getExtraPath(CONTRAST_AVERAGES_FILE), 'r') as file:
+            # Read the lines of the file
+            lines = file.readlines()
+        # Open the file for writing
+        with open(self._getExtraPath(CONTRAST_AVERAGES_FILE), "w") as file:
+            # Iterate through the lines
+            for line in lines:
+                # Replace "data_" with "data_particles" if found
+                modifiedLine = line.replace("data_", "data_particles")
+                # Write the modified line to the file
+                file.write(modifiedLine)
+
+        with open(self._getExtraPath(AVERAGES_IMAGES_FILE), 'r') as file:
+            lines = file.readlines()
+
+            # Find the index of the last non-empty line
+        lastNonEmptyInd = len(lines) - 1
+        while lastNonEmptyInd >= 0 and lines[lastNonEmptyInd].strip() == "":
+            lastNonEmptyInd -= 1
+
+        # Modify the lines
+        modifiedLines = []
+        for line in lines[:lastNonEmptyInd + 1]:
+            # Replace "data_" with "data_particles" if found
+            modifiedLine = line.replace("data_Particles", "data_particles")
+            modifiedLines.append(modifiedLine)
+
+        # Write the modified lines back to the file
+        with open(self._getExtraPath(AVERAGES_IMAGES_FILE), "w") as file:
+            file.writelines(modifiedLines)
+
+    def _loadClassesInfo(self, filename):
+        """ Read some information about the produced 2D classes
+        from the metadata file.
+        """
+        self._classesInfo = {}  # store classes info, indexed by class id
+
+        mdFileName = '%s@%s' % ('particles', filename)
+        table = emtable.Table(fileName=filename)
+
+        for classNumber, row in enumerate(table.iterRows(mdFileName)):
+            index, fn = xmippToLocation(row.get(XMIPPCOLUMNS.image.value))
+            # Store info indexed by id, we need to store the row.clone() since
+            # the same reference is used for iteration
+            self._classesInfo[classNumber + 1] = (index, fn, row)
+        self._numClass = index
+
+    def _fillClassesFromLevel(self, clsSet, update=False):
+        """ Create the SetOfClasses2D from a given iteration. """
+        self._createModelFile()
+        self._loadClassesInfo(self._getExtraPath(CONTRAST_AVERAGES_FILE))
+        # self._loadClassesInfo(self._getExtraPath('classes_classes.star'))
+        mdIter = emtable.Table.iterRows('particles@' + self._getExtraPath(AVERAGES_IMAGES_FILE))
+
+        params = {}
+        if update:
+            self.info(r'Last particle processed id is %s' % self.lastParticleProcessedId)
+            params = {"where": "id > %s" % self.lastParticleProcessedId}
+
+        # the particle with orientation parameters (all_parameters)
+        with self._lock:
+            clsSet.classifyItems(updateItemCallback=self._updateParticle,
+                                 updateClassCallback=self._updateClass,
+                                 itemDataIterator=mdIter,  # relion style
+                                 iterParams=params)
+    
+# --------------------------- Static functions --------------------------------
+def rowToAlignmentEmtable(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 alignmentRow.hasAnyColumn(ALIGNMENT_DICT.values()):
+        alignment = Transform()
+        angles = np.zeros(3)
+        shifts = np.zeros(3)
+        flip = alignmentRow.get(XMIPPCOLUMNS.flip.value, default=0.)
+
+        shifts[0] = alignmentRow.get(XMIPPCOLUMNS.shiftX.value, default=0.)
+        shifts[1] = alignmentRow.get(XMIPPCOLUMNS.shiftY.value, default=0.)
+
+        if not is2D:
+            angles[0] = alignmentRow.get(XMIPPCOLUMNS.angleRot.value, default=0.)
+            angles[1] = alignmentRow.get(XMIPPCOLUMNS.angleTilt.value, default=0.)
+            angles[2] = alignmentRow.get(XMIPPCOLUMNS.anglePsi.value, default=0.)
+            shifts[2] = alignmentRow.get(XMIPPCOLUMNS.shiftZ.value, default=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.get(XMIPPCOLUMNS.anglePsi.value, default=0.)
+            rot = alignmentRow.get(XMIPPCOLUMNS.angleRot.value, default=0.)
+            if not np.isclose(rot, 0., atol=1e-6 ) and not np.isclose(psi, 0., atol=1e-6):
+                print("HORROR rot and psi are different from zero in 2D case")
+
+            angles[0] = psi + rot
+
+        M = matrixFromGeometry(shifts, angles, inverseTransform)
+
+        if flip:
+            if alignType == ALIGN_2D:
+                M[0, :2] *= -1.  # invert only the first two columns
+                # keep x
+                M[2, 2] = -1.  # set 3D rot
+            elif alignType == ALIGN_3D:
+                M[0, :3] *= -1.  # now, invert first line excluding x
+                M[3, 3] *= -1.
+            elif alignType == ALIGN_PROJ:
+                pass
+
+        alignment.setMatrix(M)
+
+    else:
+        alignment = None
+
+    return alignment
diff --git a/xmipp3/protocols/protocol_classify_pca_streaming.py b/xmipp3/protocols/protocol_classify_pca_streaming.py
new file mode 100644
index 000000000..be0e2beb0
--- /dev/null
+++ b/xmipp3/protocols/protocol_classify_pca_streaming.py
@@ -0,0 +1,494 @@
+# ******************************************************************************
+# *
+# * Authors: Erney Ramirez Aportela (eramirez@cnb.csic.es)
+# *          Daniel Marchan Torres (da.marchan@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'
+# *
+# ******************************************************************************
+import sys
+import emtable
+import os
+from datetime import datetime
+import time
+
+from pyworkflow.utils import prettyTime
+from pyworkflow import VERSION_3_0
+from pyworkflow.object import Set
+from pyworkflow.protocol.params import IntParam
+from pyworkflow.protocol import ProtStreamingBase, STEPS_PARALLEL
+from pyworkflow.constants import BETA
+
+from pwem.objects import SetOfClasses2D, SetOfAverages
+from pwem.constants import ALIGN_NONE, ALIGN_2D
+
+from xmipp3.convert import (readSetOfParticles, writeSetOfParticles,
+                            writeSetOfClasses2D)
+from xmipp3.protocols.protocol_classify_pca import updateEnviron, XmippProtClassifyPca
+
+OUTPUT_CLASSES = "outputClasses"
+OUTPUT_AVERAGES = "outputAverages"
+PCA_FILE = "pca_done.txt"
+CLASSIFICATION_FILE = "classification_done.txt"
+LAST_DONE_FILE = "last_done.txt"
+BATCH_UPDATE = 5000
+
+
+class XmippProtClassifyPcaStreaming(ProtStreamingBase, XmippProtClassifyPca):
+    """ Classifies a set of images. """
+
+    _label = '2D classification pca streaming'
+    _lastUpdateVersion = VERSION_3_0
+    _conda_env = 'xmipp_pyTorch'
+    _devStatus = BETA
+
+    # Mode
+    CREATE_CLASSES = 0
+    UPDATE_CLASSES = 1
+
+    _possibleOutputs = {OUTPUT_CLASSES: SetOfClasses2D,
+                        OUTPUT_AVERAGES: SetOfAverages}
+
+    def __init__(self, **args):
+        XmippProtClassifyPca.__init__(self, **args)
+        self.stepsExecutionMode = STEPS_PARALLEL
+
+    # --------------------------- DEFINE param functions ------------------------
+    def _defineParams(self, form):
+        form = self._defineCommonParams(form)
+        form.addSection(label='Classification')
+        form.addParam('classificationBatch', IntParam, default=50000,
+                      condition="not mode",
+                      label="particles for initial classification",
+                      help='Number of particles for an initial classification to compute the 2D references')
+
+        form.addParallelSection(threads=3, mpi=1)
+
+    # --------------------------- INSERT steps functions ----------------------
+    def stepsGeneratorStep(self) -> None:
+        """
+        This step should be implemented by any streaming protocol.
+        It should check its input and when ready conditions are met
+        call the self._insertFunctionStep method.
+        """
+        self._initialStep()
+        self.newDeps = []
+        newParticlesSet = self._loadEmptyParticleSet()
+
+        if self.isContinued() and self._isPcaDone():
+            self.info('Continue protocol')
+            self._updateVarsToContinue()
+
+        while not self.finish:
+            if not self._newParticlesToProcess():
+                self.info('No new particles')
+            else:
+                with self._lock:  # This lock needs to be here since the classifyItems access the inputSet
+                    # and can create a race condition
+                    particlesSet = self._loadInputParticleSet()
+                    self.streamState = particlesSet.getStreamState()
+                    where = None
+                    if self.lastCreationTime:
+                        where = 'creation>"' + str(self.lastCreationTime) + '"'
+
+                    for particle in particlesSet.iterItems(orderBy='creation', direction='ASC', where=where):
+                        tmp = particle.getObjCreation()
+                        newParticlesSet.append(particle.clone())
+
+                self.lastCreationTime = tmp
+                self.info('%d new particles' % len(newParticlesSet))
+                self.info('Last creation time REGISTER %s' % self.lastCreationTime)
+
+                # ------------------------------------ PCA TRAINING ----------------------------------------------
+                if self._doPcaTraining(newParticlesSet):
+                    self.pcaStep = self._insertFunctionStep(self.runPCASteps, newParticlesSet, prerequisites=[])
+                    if len(newParticlesSet) == len(self.inputParticles.get()) and self.streamState == Set.STREAM_CLOSED:
+                        self.staticRun = True
+                        self.info('Static Run')
+                # ------------------------------------ CLASSIFICATION ----------------------------------------------
+                if self._doClassification(newParticlesSet):
+                    self._insertClassificationSteps(newParticlesSet, self.lastCreationTime)
+                    newParticlesSet = self._loadEmptyParticleSet()
+
+            if self.streamState == Set.STREAM_CLOSED:
+                self.info('Stream closed')
+                # Finish everything and close output sets
+                if len(newParticlesSet):
+                    self.info('Finish processing with last batch %d' % len(newParticlesSet))
+                    self.lastRound = True
+                else:
+                    self._insertFunctionStep(self.closeOutputStep, prerequisites=self.newDeps)
+                    self.finish = True
+                    continue  # To avoid waiting 1 min
+
+            sys.stdout.flush()
+            time.sleep(60)
+
+        sys.stdout.flush()  # One last flush
+
+    # --------------------------- STEPS functions -------------------------------
+    def _initialStep(self):
+        self.finish = False
+        self.lastCreationTime = 0
+        self.lastCreationTimeProcessed = 0
+        self.streamState = Set.STREAM_OPEN
+        self.lastRound = False
+        self.pcaLaunch = False
+        self.classificationLaunch = False
+        self.classificationRound = 0
+        self.firstTimeDone = False
+        self.staticRun = False
+        # Initialize files
+        self._initFnStep()
+
+    def _initFnStep(self):
+        updateEnviron(self.gpuList.get())
+        self.imgsPcaXmd = self._getExtraPath('images_pca.xmd')
+        self.imgsPcaXmdOut = self._getTmpPath('images_pca.xmd')  # Wiener
+        self.imgsPcaFn = self._getTmpPath('images_pca.mrc')
+        self.imgsOrigXmd = self._getExtraPath('imagesInput_.xmd')
+        self.imgsXmd = self._getTmpPath('images_.xmd')
+        self.imgsFn = self._getTmpPath('images_.mrc')
+        self.refXmd = self._getTmpPath('references.xmd')
+        self.ref = self._getExtraPath('classes.mrcs')
+        self.sigmaProt = self.sigma.get()
+        if self.sigmaProt == -1:
+            self.sigmaProt = self.inputParticles.get().getDimensions()[0] / 3
+
+        self.sampling = self.inputParticles.get().getSamplingRate()
+        resolution = self.resolution.get()
+        if resolution < 2 * self.sampling:
+            resolution = (2 * self.sampling) + 0.5
+        self.resolutionPca = resolution
+
+        if self.mode == self.UPDATE_CLASSES:
+            self.numberClasses = len(self.initialClasses.get())
+            self.classificationBatch.set(BATCH_UPDATE)
+        else:
+            self.numberClasses = self.numberOfClasses.get()
+
+    def runPCASteps(self, newParticlesSet):
+        # Run PCA steps
+        self.pcaLaunch = True
+        if self.correctCtf:
+            self.convertInputStep(newParticlesSet, self.imgsPcaXmd, self.imgsPcaXmdOut)  # Wiener filter
+        else:
+            self.convertInputStep(newParticlesSet, self.imgsPcaXmd, self.imgsPcaFn)
+        numTrain = min(len(newParticlesSet), self.training.get())
+        self.pcaTraining(self.imgsPcaFn, self.resolutionPca, numTrain)
+        self.pcaLaunch = False
+        self._setPcaDone()
+
+    def _insertClassificationSteps(self, newParticlesSet, lastCreationTime):
+        self._updateFnClassification()
+        classStep = self._insertFunctionStep(self.runClassificationSteps,
+                                             newParticlesSet, prerequisites=self.pcaStep)
+        updateStep = self._insertFunctionStep(self.updateOutputSetOfClasses,
+                                              lastCreationTime, Set.STREAM_OPEN, prerequisites=classStep)
+        self.newDeps.append(updateStep)
+
+    def runClassificationSteps(self, newParticlesSet):
+        self.classificationLaunch = True
+        self.convertInputStep(newParticlesSet, self.imgsOrigXmd, self.imgsFn)
+        self.classification(self.imgsFn, self.numberClasses,
+                            self.imgsOrigXmd, self.mask.get(), self.sigmaProt)
+        self.classificationLaunch = False
+
+    def convertInputStep(self, input, outputOrig, outputMRC):
+        if self._isPcaDone() and self.staticRun:  # Static run: reuse files to avoid doing this twice
+            self.imgsOrigXmd = self.imgsPcaXmd
+            self.imgsFn = self.imgsPcaFn
+        else:
+            writeSetOfParticles(input, outputOrig)
+
+            if self.correctCtf:
+                # args = ' -i  %s -o %s --sampling_rate %s ' % (outputOrig, outputMRC, self.sampling)
+                # self.runJob("xmipp_ctf_correct_wiener2d", args,  numberOfMpi=8)
+                # ------------ WIENER -----------------------
+                args = (' -i %s -o %s --pixel_size %s --spherical_aberration %s '
+                       '--voltage %s --q0 %s --batch 512 --padding 2 --device cuda:%d') % \
+                      (outputOrig, outputMRC, self.sampling, self.acquisition.getSphericalAberration(),
+                       self.acquisition.getVoltage(), self.acquisition.getAmplitudeContrast(), int(self.gpuList.get()))
+                env = self.getCondaEnv()
+                env = self._setEnvVariables(env)
+                self.runJob("xmipp_swiftalign_wiener_2d", args, numberOfMpi=1, env=env)
+            else:
+                args = ' -i  %s -o %s  ' % (outputOrig, outputMRC)
+                self.runJob("xmipp_image_convert", args, numberOfMpi=1)
+
+        # For classification update
+        if self.mode == self.UPDATE_CLASSES and not self.firstTimeDone:
+            initial2DClasses = self.initialClasses.get()
+            self.info('Write classes')
+            if isinstance(initial2DClasses, SetOfClasses2D):
+                writeSetOfClasses2D(initial2DClasses,
+                                    self.refXmd, writeParticles=False)
+            else:
+                writeSetOfParticles(initial2DClasses,
+                                    self.refXmd,
+                                    alignType=ALIGN_NONE)
+
+            args = ' -i  %s -o %s  ' % (self.refXmd, self.ref)
+            self.runJob("xmipp_image_convert", args, numberOfMpi=1)
+            self.firstTimeDone = True
+
+    def pcaTraining(self, inputIm, resolutionTrain, numTrain):
+        args = ' -i %s  -s %s -hr %s -lr 530 -p %s -t %s -o %s/train_pca  --batchPCA' % \
+               (inputIm, self.sampling, resolutionTrain, self.coef.get(), numTrain, self._getExtraPath())
+
+        env = self.getCondaEnv()
+        env = self._setEnvVariables(env)
+        self.runJob("xmipp_classify_pca_train", args, numberOfMpi=1, env=env)
+
+    def classification(self, inputIm, numClass, stfile, mask, sigma):
+        args = ' -i %s -c %s -b %s/train_pca_bands.pt -v %s/train_pca_vecs.pt -o %s/classes -stExp %s' % \
+               (inputIm, numClass, self._getExtraPath(), self._getExtraPath(), self._getExtraPath(),
+                stfile)
+        if mask:
+            args += ' --mask --sigma %s ' % (sigma)
+
+        if self.mode == self.UPDATE_CLASSES or self._isClassificationDone():
+            args += ' -r %s ' % self.ref
+
+        env = self.getCondaEnv()
+        env = self._setEnvVariables(env)
+        self.runJob("xmipp_classify_pca", args, numberOfMpi=1, env=env)
+
+    def updateOutputSetOfClasses(self, lastCreationTime, streamMode):
+        outputName = OUTPUT_CLASSES
+        outputClasses, update = self._loadOutputSet(outputName)
+
+        self._fillClassesFromLevel(outputClasses, update)
+        self._updateOutputSet(outputName, outputClasses, streamMode)
+        self._updateOutputAverages(update)
+
+        if not update:  # First time
+            self._defineSourceRelation(self._getInputPointer(), outputClasses)
+            self._setClassificationDone()
+            self.numberClasses = len(outputClasses)  # In case the original number of classes is not reached
+
+        self.lastCreationTimeProcessed = lastCreationTime
+        self.info(r'Last creation time processed UPDATED is %s' % str(self.lastCreationTimeProcessed))
+        self._writeLastDone(str(self.lastCreationTimeProcessed))
+        self._writeLastClassificationRound(self.classificationRound)
+
+        self.info(r'Last classification round processed is %d' % self.classificationRound)
+        self.classificationRound += 1
+
+    def _updateOutputAverages(self, update):
+        outRefs = self._loadOutputAverageSet()
+        readSetOfParticles(self.ref, outRefs)
+        self._updateOutputSet(OUTPUT_AVERAGES, outRefs, Set.STREAM_CLOSED)
+        if not update:  # First Time
+            self._defineSourceRelation(self._getInputPointer(), outRefs)
+
+    def closeOutputStep(self):
+        self._closeOutputSet()
+
+    # --------------------------- UTILS functions -----------------------------
+    def _loadInputParticleSet(self):
+        """ Returns te input set of particles"""
+        partSet = self.inputParticles.get()
+        partSet.loadAllProperties()
+
+        return partSet
+
+    def _getInputPointer(self):
+        return self.inputParticles
+
+    def _loadEmptyParticleSet(self):
+        partSet = self.inputParticles.get()
+        partSet.loadAllProperties()
+        self.acquisition = partSet.getAcquisition()
+        copyPartSet = self._createSetOfParticles()
+        copyPartSet.copyInfo(partSet)
+        copyPartSet.setSamplingRate(self.sampling)
+        copyPartSet.setAcquisition(self.acquisition)
+
+        return copyPartSet
+
+    def _setEnvVariables(self, env):
+        """ Method to set all the environment variables needed to run PCA program """
+        env['LD_LIBRARY_PATH'] = ''
+        # Limit the number of threads
+        env['OMP_NUM_THREADS'] = '12'
+        env['MKL_NUM_THREADS'] = '12'
+        return env
+
+    def _updateFnClassification(self):
+        """Update input based on the iteration it is"""
+        self.imgsOrigXmd = updateFileName(self.imgsOrigXmd, self.classificationRound)
+        self.imgsXmd = updateFileName(self.imgsXmd, self.classificationRound)
+        self.imgsFn = updateFileName(self.imgsFn, self.classificationRound)
+        self.info('Starts classification round: %d' % self.classificationRound)
+
+    def _newParticlesToProcess(self):
+        particlesFile = self.inputParticles.get().getFileName()
+        now = datetime.now()
+        self.lastCheck = getattr(self, 'lastCheck', now)
+        mTime = datetime.fromtimestamp(os.path.getmtime(particlesFile))
+        self.debug('Last check: %s, modification: %s'
+                   % (self.lastCheck,
+                      prettyTime(mTime)))
+        # If the input have not changed since our last check,
+        # it does not make sense to check for new input data
+        if (self.lastCheck > mTime and self.lastCreationTime) and not self.lastRound:
+            newParticlesBool = False
+        else:
+            newParticlesBool = True
+
+        self.lastCheck = now
+        return newParticlesBool
+
+    def _fillClassesFromLevel(self, clsSet, update=False):
+        """ Create the SetOfClasses2D from a given iteration. """
+        self._createModelFile()
+        self._loadClassesInfo(self._getExtraPath('classes_contrast_classes.star'))
+        mdIter = emtable.Table.iterRows('particles@' + self._getExtraPath('classes_images.star'))
+
+        params = {}
+        if update:
+            self.info(r'Last creation time processed is %s' % str(self.lastCreationTimeProcessed))
+            params = {"where": 'creation>"' + str(self.lastCreationTimeProcessed) + '"'}
+
+        with self._lock:
+            clsSet.classifyItems(updateItemCallback=self._updateParticle,
+                                 updateClassCallback=self._updateClass,
+                                 itemDataIterator=mdIter,  # relion style
+                                 iterParams=params,
+                                 doClone=False,  # So the creation time is maintained
+                                 raiseOnNextFailure=False)  # So streaming can happen
+
+    def _loadOutputSet(self, outputName):
+        """
+        Load the output set if it exists or create a new one.
+        """
+        outputSet = getattr(self, outputName, None)
+        update = False
+        if outputSet is None:
+            outputSet = self._createSetOfClasses2D(self._getInputPointer())
+            outputSet.setStreamState(Set.STREAM_OPEN)
+        else:
+            update = True
+
+        return outputSet, update
+
+    def _loadOutputAverageSet(self):
+        """
+        Load an empty output setOfAverages
+        """
+        outputRefs = self._createSetOfAverages()  # We need to create always an empty set since we need to rebuild it
+        outputRefs.copyInfo(self.inputParticles.get())
+        outputRefs.setSamplingRate(self.sampling)
+        outputRefs.setAlignment(ALIGN_2D)
+
+        return outputRefs
+
+    def _doPcaTraining(self, newParticlesSet):
+        """ Two cases for launching PCA steps:
+            - If there are enough particles and PCA has not been launched or finished.
+            - Launch if it's the last round of new particles and PCA has not been launched or finished. """
+        return (len(newParticlesSet) >= self.training.get() and not self._isPcaDone() and not self.pcaLaunch) or \
+            (self.lastRound and not self._isPcaDone() and not self.pcaLaunch)
+
+    def _doClassification(self, newParticlesSet):
+        """ Three cases for launching Classification
+            - First round of classification: enough particles and PCA done
+            - Update classification: classification done, PCA done and enough batch size
+            - First or update classification with a smaller batch: last round of particles and not classification launch
+        """
+        return (len(newParticlesSet) >= self.classificationBatch.get() and self._isPcaDone()
+                and not self._isClassificationDone() and not self.classificationLaunch) \
+            or (len(newParticlesSet) >= BATCH_UPDATE and self._isClassificationDone() and self._isPcaDone()
+                and not self.classificationLaunch) \
+            or (self.lastRound and not self.classificationLaunch)
+
+    def _isPcaDone(self):
+        done = False
+        if os.path.exists(self._getExtraPath(PCA_FILE)):
+            done = True
+        return done
+
+    def _setPcaDone(self):
+        with open(self._getExtraPath(PCA_FILE), "w") as file:
+            file.write('%d' % self.pcaStep)
+
+    def _getPcaStep(self):
+        " If continue then use this to collect the code for the PCA step"
+        with open(self._getExtraPath(PCA_FILE), "r") as file:
+            content = file.read()
+            return int(content)
+
+    def _isClassificationDone(self):
+        done = False
+        if os.path.exists(self._getExtraPath(CLASSIFICATION_FILE)):
+            done = True
+
+        if self.mode == self.UPDATE_CLASSES:
+            done = True
+
+        return done
+
+    def _setClassificationDone(self):
+        with open(self._getExtraPath(CLASSIFICATION_FILE), "w"):
+            self.debug("Creating Classification DONE file")
+
+    def _writeLastClassificationRound(self, classificationRound):
+        with open(self._getExtraPath(CLASSIFICATION_FILE), "w") as file:
+            file.write('%d' % classificationRound)
+
+    def _getLastClassificationRound(self):
+        with open(self._getExtraPath(CLASSIFICATION_FILE), "r") as file:
+            content = file.read()
+            return int(content)
+
+    def _writeLastDone(self, creationTime):
+        """ Write to a text file the last item creation time done. """
+        with open(self._getExtraPath(LAST_DONE_FILE), 'w') as file:
+            file.write('%s' % creationTime)
+
+    def _getLastDone(self):
+        # Open the file in read mode and read the number
+        with open(self._getExtraPath(LAST_DONE_FILE), "r") as file:
+            content = file.read()
+        return str(content)
+
+    def _updateVarsToContinue(self):
+        """ Method to if needed and the protocol is set to continue then it will see in which state it was stopped """
+
+        self.pcaStep = []
+
+        if self._isClassificationDone():
+            self.lastCreationTime = self._getLastDone()
+            self.classificationRound = self._getLastClassificationRound() + 1  # Since this is the last processed
+        else:
+            self.lastCreationTime = 0
+            self.classificationRound = 0
+
+        self.lastCreationTimeProcessed = self.lastCreationTime
+        # Convert the string to a datetime object
+        self.lastCheck = datetime.strptime(self.lastCreationTime, '%Y-%m-%d %H:%M:%S')
+
+# --------------------------- Static functions --------------------------------
+def updateFileName(filepath, round):
+    filename = os.path.basename(filepath)
+    newFilename = f"{filename[:filename.find('_')]}_{round}{filename[filename.rfind('.'):]}"
+    return os.path.join(os.path.dirname(filepath), newFilename)
diff --git a/xmipp3/protocols/protocol_compare_reprojections.py b/xmipp3/protocols/protocol_compare_reprojections.py
index 03eb274c9..c2c297d17 100644
--- a/xmipp3/protocols/protocol_compare_reprojections.py
+++ b/xmipp3/protocols/protocol_compare_reprojections.py
@@ -3,6 +3,7 @@
 # * Authors:     Josue Gomez Blanco (josue.gomez-blanco@mcgill.ca)
 # *              Estrella Fernandez Gimenez (me.fernandez@cnb.csic.es)
 # *              (produce residuals improved with projection subtraction)
+# *              Daniel Marchan Torres (da.marchan@cnb.csic.es)
 # *
 # * Unidad de  Bioinformatica of Centro Nacional de Biotecnologia , CSIC
 # *
@@ -25,21 +26,30 @@
 # *  e-mail address 'scipion@cnb.csic.es'
 # *
 # **************************************************************************
-
 from math import floor
 import os
+import json
+import numpy as np
 
-from pyworkflow import VERSION_1_1
-from pyworkflow.protocol.params import (PointerParam, StringParam, FloatParam, BooleanParam, IntParam)
-from pyworkflow.utils.path import cleanPath
 from pwem.protocols import ProtAnalysis3D
-from pwem.objects import (SetOfClasses2D, Image, SetOfAverages, SetOfParticles)
+from pwem.objects import (SetOfClasses2D, Image, SetOfAverages, SetOfParticles,
+                          SetOfVolumes, SetOfClasses3D, Volume, EMSet, Class3D)
 import pwem.emlib.metadata as md
+from pwem.emlib.image import ImageHandler
+from pwem import emlib
+
 from pyworkflow.protocol.constants import LEVEL_ADVANCED
+from pyworkflow import UPDATED
+from pyworkflow import VERSION_3_0
+from pyworkflow.protocol.params import (PointerParam, StringParam, FloatParam, BooleanParam, EnumParam)
+from pyworkflow.utils.path import cleanPath, cleanPattern
 
-from pwem import emlib
 from xmipp3.base import ProjMatcher
-from xmipp3.convert import setXmippAttributes, xmippToLocation
+from xmipp3.convert import setXmippAttributes, xmippToLocation, getXmippAttribute
+from ..convert import writeSetOfClasses2D, writeSetOfParticles
+
+NEW_SAMPLING_RATE = 3.0
+OUTPUTS_FN = "outputs.txt"
 
         
 class XmippProtCompareReprojections(ProtAnalysis3D, ProjMatcher):
@@ -52,8 +62,15 @@ class XmippProtCompareReprojections(ProtAnalysis3D, ProjMatcher):
     values particularly low or high, may indicate outliers."""
 
     _label = 'compare reprojections'
-    _lastUpdateVersion = VERSION_1_1
+    _lastUpdateVersion = VERSION_3_0
+    _devStatus = UPDATED
+    _possibleOutputs = {}
+
+    PARTICLES = 0
+    VOLUME = 1
+    BOTH = 2
     
+
     def __init__(self, **args):
         ProtAnalysis3D.__init__(self, **args)
         self._classesInfo = dict()
@@ -61,15 +78,19 @@ def __init__(self, **args):
     # --------------------------- DEFINE param functions --------------------------------------------
     def _defineParams(self, form):
         form.addSection(label='Input')
-        form.addParam('inputSet', PointerParam, label="Input images", important=True, 
-                      pointerClass='SetOfClasses2D, SetOfAverages, SetOfParticles')
-        form.addParam('inputVolume', PointerParam, label="Volume to compare images to", important=True,
-                      pointerClass='Volume', help='Volume to be used for class comparison')
+        form.addParam('inputSet2D', PointerParam, label="Input images", important=True,
+                      pointerClass='SetOfClasses2D, SetOfClasses3D, SetOfAverages, SetOfParticles')
+        form.addParam('inputSet3D', PointerParam, label="Volume to compare images to", important=True,
+                      pointerClass='Volume, SetOfVolumes, SetOfClasses3D',
+                      help='Volume to be used for class comparison')
         form.addParam('useAssignment', BooleanParam, default=True, label='Use input angular assignment (if available)')
         form.addParam('optimizeGray', BooleanParam, default=False, label='Optimize gray scale')
         form.addParam('ignoreCTF', BooleanParam, default=True, label='Ignore CTF',
                       help='By ignoring the CTF you will create projections more similar to what a person expects, '
                            'while by using the CTF you will create projections more similar to what the microscope sees')
+        form.addParam('doDownSample', BooleanParam, default=True, label='Downsample',
+                      help='If accepted the input volumes and the input 2d classes will be downsample to 3 A/px. '
+                           'This will help to reduce the computation time.')
         form.addParam('doEvaluateResiduals', BooleanParam, default=False, expertLevel=LEVEL_ADVANCED,
                       label='Evaluate residuals',
                       help='If this option is chosen, then the residual covariance matrix is calculated and '
@@ -87,111 +108,275 @@ def _defineParams(self, form):
         form.addParam('sigma', FloatParam, label="Decay of the filter (sigma): ", default=3, condition='resol',
                       help='Decay of the filter (sigma parameter) to smooth the mask transition',
                       expertLevel=LEVEL_ADVANCED)
+        form.addSection(label='Output')
+        form.addParam('doRanking', BooleanParam, default=True, label='Rank the set of Volumes/3D Classes',
+                      help='If accepted please check that you have an input set with several 3D objects '
+                           '(Classes 3D or Volumes). It will rank all the Volume 3D Class/Volume.')
+        form.addParam('doExtraction', BooleanParam, default=True, label='Extract best Volume/3D Class',
+                      condition='doRanking',
+                      help='If accepted please check that you have an input set with several 3D objects '
+                           '(Classes 3D or Volumes). It will extract the Volume and/or the Particles from '
+                           'the best rated 3D Class/Volume.')
+        form.addParam('extractOption', EnumParam,
+                      choices=['Particles', 'Volume', 'Both'],
+                      default=self.BOTH,
+                      condition='doExtraction and doRanking',
+                      label="Extraction option", display=EnumParam.DISPLAY_COMBO,
+                      help='Select an option to extract from the 3D Classes: \n '
+                           '_Particles_: Extract the set of particles from the selected class. \n '
+                           '_Volume_: Extract the volume from the selected class. \n'
+                           '_Both_: Extract the volume and particles from the selected class.')
+
         form.addParallelSection(threads=0, mpi=8)
     
     # --------------------------- INSERT steps functions --------------------------------------------
     def _insertAllSteps(self):
         """ Convert input images if necessary """
-        self.imgsFn = self._getExtraPath('residuals.xmd')
-        vol = self.inputVolume.get()
-        self._insertFunctionStep("convertStep", self.imgsFn)
+        self.initialStep()
+        self._insertFunctionStep(self.convertStep)
+
+        if self.doDownSample:
+            self._insertFunctionStep(self.downSampleStep, self.imgsOrigFn)
 
         if self._useProjMatching():
-            anglesFn = self._getExtraPath('angles.xmd')
-            self._insertFunctionStep("projMatchStep", self.inputVolume.get().getFileName(), self.angularSampling.get(),
-                                     self.symmetryGroup.get(), self.imgsFn, anglesFn, self.inputVolume.get().getDim()[0])
+            self._insertFunctionStep(self.insertProjMatchStep, self.fnVolDict, self.angularSampling.get(),
+                                     self.symmetryGroup.get(), self.anglesDict, self.galleryDict)
         else:
-            anglesFn = self.imgsFn
+            tmpDict = self.anglesDict
+            self.anglesDict = {volId: self.imgsFn for volId in tmpDict}
 
-        self._insertFunctionStep("produceResiduals", vol.getFileName(), anglesFn, vol.getSamplingRate())
+        self._insertFunctionStep(self.produceResiduals, self.anglesDict,
+                                 NEW_SAMPLING_RATE if self.doDownSample else self.inputSet3D.get().getSamplingRate())
 
         if self.doEvaluateResiduals.get():
-            self._insertFunctionStep("evaluateResiduals")
+            self._insertFunctionStep(self.evaluateResiduals)
 
-        self._insertFunctionStep("createOutputStep")
+        self._insertFunctionStep(self.createOutputStep)
 
     # --------------------------- STEPS functions ---------------------------------------------------
-    def convertStep(self, imgsFn):
-        from ..convert import writeSetOfClasses2D, writeSetOfParticles
-        imgSet = self.inputSet.get()
+    def initialStep(self):
+        self.samplingRateVol = self.inputSet3D.get().getSamplingRate()
+        self.samplingRateAverages = self.inputSet2D.get().getSamplingRate()
+        self.imgsOrigFn = self._getExtraPath('residuals.xmd')
+        self.fnVolDict = {}
+        self.anglesDict = {}
+        self.galleryDict = {}
+        self.outputNamesDict = {}
+        volSet = self.inputSet3D.get()
+        self.xDim = self._getWorkingDimensions(self.samplingRateAverages, self.samplingRateVol, self.doDownSample.get())
+
+        if self.doDownSample and self.samplingRateAverages < NEW_SAMPLING_RATE:
+            self.imgsFn = self._getExtraPath('residuals_downsample.xmd')
+        else:
+            self.imgsFn = self.imgsOrigFn
+
+        if isinstance(volSet, EMSet):
+            for vol in volSet.iterItems(orderBy='id', direction='ASC'):
+                vid = vol.getObjId()
+                self.fnVolDict[vid] = self._getTmpPath("volume_%d.vol" % vid)
+                self.anglesDict[vid] = self._getExtraPath('angles_%d.xmd' % vid)
+                self.galleryDict[vid] = self._getExtraPath('gallery_%d.stk' % vid)
+                self.outputNamesDict[vid] = "reprojections_vol%d" % vid
+        else:
+            vid = volSet.getObjId()
+            self.fnVolDict[vid] = self._getTmpPath("volume_%d.vol" % vid)
+            self.anglesDict[vid] = self._getExtraPath('angles_%d.xmd' % vid)
+            self.galleryDict[vid] = self._getExtraPath('gallery_%d.stk' % vid)
+            self.outputNamesDict[vid] = "reprojections_vol%d" % vid
+
+        with open(self._getExtraPath(OUTPUTS_FN), 'w'):
+            pass  # Creates output file for visualizations warnings
+
+
+    def convertStep(self):
+        # Convert input images (SetOfClasses2D, SetOfAverages and SetOfParticles)
+        imgSet = self.inputSet2D.get()
+        imgsFn = self.imgsOrigFn
         if isinstance(imgSet, SetOfClasses2D):
-            writeSetOfClasses2D(imgSet, self.imgsFn, writeParticles=True)
+            writeSetOfClasses2D(imgSet, imgsFn, writeParticles=True)
         else:
-            writeSetOfParticles(imgSet, self.imgsFn)
-        from pwem.emlib.image import ImageHandler
-        img = ImageHandler()
-        fnVol = self._getTmpPath("volume.vol")
-        img.convert(self.inputVolume.get(), fnVol)
-        xdim = self.inputVolume.get().getDim()[0]
-        if xdim != self._getDimensions():
-            self.runJob("xmipp_image_resize", "-i %s --dim %d" % (fnVol, self._getDimensions()), numberOfMpi=1)
-
-    def produceResiduals(self, fnVol, fnAngles, Ts):
-        anglesOutFn = self._getExtraPath("anglesCont.stk")
-        projectionsOutFn = self._getExtraPath("projections.stk")
-        fnVol = self._getTmpPath("volume.vol")
-        xdim = self._getDimensions()
-        args = "-i %s -o %s --ref %s --optimizeAngles --optimizeShift --max_shift %d --oprojections %s --sampling %f" \
-               % (fnAngles, anglesOutFn, fnVol, floor(xdim*0.05), projectionsOutFn, Ts)
-        self.fnResiduals = self._getExtraPath("residuals.stk")
-        if self.doEvaluateResiduals:
-            args += " --oresiduals %s" % self.fnResiduals
+            writeSetOfParticles(imgSet, imgsFn)
+        # Convert input volumes (SetOfClasses3D, SetOfVolumes and Volume)
+        volSet = self.inputSet3D.get()
+        for volId, fnVol in self.fnVolDict.items():
+            if isinstance(volSet, Volume):
+                vol = volSet
+            else:
+                volCl = volSet.getItem("id", volId)
+                if isinstance(volCl, Class3D):
+                    vol = volCl.getRepresentative().clone()
+                else:
+                    vol = volCl
+            img = ImageHandler()
+            fnVol = self.fnVolDict[volId]
+            self.info("Registering volume %s" % fnVol)
+            img.convert(vol, fnVol)
 
-        if self.ignoreCTF:
-            args += " --ignoreCTF "
-        if self.optimizeGray:
-            args += "--optimizeGray --max_gray_scale 0.95 "
-        self.runJob("xmipp_angular_continuous_assign2", args)
+            # In case input volume does not match 2D references size (If downsample is activated this is not needed)
+            xdimVol = vol.getDim()[0]
+            xdimImg = self._getDimensionsImages()
+            if xdimVol != xdimImg and not self.doDownSample:
+                self.xDim = xdimImg
+                self.runJob("xmipp_image_resize", "-i %s --dim %d"
+                            % (fnVol, self.xDim), numberOfMpi=1)
 
-        if self.doEvaluateResiduals:
-            self._computeResiduals(fnVol)
+    def downSampleStep(self, imgsOrigFn):
+        # Calculate new sampling rate
+        if NEW_SAMPLING_RATE < self.samplingRateVol or NEW_SAMPLING_RATE < self.samplingRateAverages:
+            newSamplingRate = max(self.samplingRateVol, self.samplingRateAverages)
+            self.info("The target sampling rate is smaller than the inputs sampling rate, new target sampling rate "
+                      "is set to the biggest one of both inputs %f" % newSamplingRate)
+        else:
+            newSamplingRate = NEW_SAMPLING_RATE
+
+        # Downsample the average 2D classes
+        if self.samplingRateAverages < newSamplingRate:
+            imgsFn = self._getExtraPath('residuals_downsample.mrcs')
+            self.info('Downsampling the input images...')
+            self.runJob("xmipp_image_resize", "-i %s -o %s --dim %d"
+                        % (imgsOrigFn, imgsFn, self.xDim), numberOfMpi=1)
+        else:
+            self.info('The target sampling rate %f <= %f the current one, skipping downsample for the 2D Images.'
+                      % (newSamplingRate, self.samplingRateAverages))
+
+        # Downsample Volume/s
+        if self.samplingRateVol < newSamplingRate:
+            self.info('Downsampling the input volumes...')
+            for fnVol in self.fnVolDict.values():
+                self.runJob("xmipp_image_resize", "-i %s --dim %d"  # Should we store somewhere the resized Vol?
+                            % (fnVol, self.xDim), numberOfMpi=self.numberOfMpi)
+        else:
+            self.info('The target sampling rate %f <= %f the current one, skipping downsample for the 3D Volumes.'
+                      % (newSamplingRate, self.samplingRateVol))
+
+    def insertProjMatchStep(self, volumesDict, angularSampling, symmetryGroup, fnAnglesDict, galleryDict):
+        xDim = self.xDim
+        images = self.imgsFn
+        for volId, fnVol in volumesDict.items():
+            self.info('Generating a gallery of projections for volume %s' % fnVol)
+            fnAngles = fnAnglesDict[volId]
+            fnGallery = galleryDict[volId]
+            self.projMatchStepAdapted(fnVol, angularSampling, symmetryGroup, images, fnAngles, fnGallery, xDim, volId)
+
+    def projMatchStepAdapted(self, volume, angularSampling, symmetryGroup, images, fnAngles, fnGallery, xDim, volId):
+        # Generate gallery of projections
+        if volume.endswith('.mrc'):
+            volume += ":mrc"
+
+        self.runJob("xmipp_angular_project_library",
+                    "-i %s -o %s --sampling_rate %f --sym %s --method fourier 1 0.25 bspline "
+                    "--compute_neighbors --angular_distance -1 --experimental_images %s"
+                    % (volume, fnGallery, angularSampling, symmetryGroup, images))
+
+        # Assign angles
+        self.runJob("xmipp_angular_projection_matching",
+                    "-i %s -o %s --ref %s --Ri 0 --Ro %s --max_shift 1000 "
+                    "--search5d_shift %s --search5d_step  %s --append"
+                    % (images, fnAngles, fnGallery, str(xDim / 2),
+                       str(int(xDim / 10)), str(int(xDim / 25))))
+
+        cleanPath(self._getExtraPath('gallery_%d_sampling.xmd' % volId))
+        cleanPath(self._getExtraPath('gallery_%d_angles.doc' % volId))
+        cleanPath(self._getExtraPath('gallery_%d.doc' % volId))
+
+        # Write angles in the original file and sort
+        MD = emlib.MetaData(fnAngles)
+        for id in MD:
+            galleryReference = MD.getValue(emlib.MDL_REF, id)
+            MD.setValue(emlib.MDL_IMAGE_REF, "%05d@%s" % (galleryReference + 1, fnGallery), id)
+        MD.write(fnAngles)
+
+    def produceResiduals(self, anglesDict, tS):
+        self.fnResidualsDict = {}
+        for volId, fnAngles in anglesDict.items():
+            anglesOutFn = self._getExtraPath("anglesCont_%d.stk" % volId)
+            projectionsOutFn = self._getExtraPath("projections_%d.stk" % volId)
+            fnVol = self.fnVolDict[volId]
+            xDim = self.xDim
+            args = "-i %s -o %s --ref %s --optimizeAngles --optimizeShift --max_shift %d --oprojections %s --sampling %f" \
+                   % (fnAngles, anglesOutFn, fnVol, floor(xDim*0.05), projectionsOutFn, tS)
+
+            fnResiduals = self._getExtraPath("residuals_%d.mrcs" % volId)
+            self.fnResidualsDict[volId] = fnResiduals
+
+            if self.doEvaluateResiduals:
+                args += " --oresiduals %s" % fnResiduals
+
+            if self.ignoreCTF:
+                args += " --ignoreCTF "
+            if self.optimizeGray:
+                args += "--optimizeGray --max_gray_scale 0.95 "
+
+            self.runJob("xmipp_angular_continuous_assign2", args)
 
     def evaluateResiduals(self):
-        # Evaluate each image
-        fnAutoCorrelations = self._getExtraPath("autocorrelations.xmd")
-        stkAutoCorrelations = self._getExtraPath("autocorrelations.stk")
-        stkResiduals = self.fnResiduals
-        anglesOutFn = self._getExtraPath("anglesCont.xmd")
-        self.runJob("xmipp_image_residuals", " -i %s -o %s --save_metadata_stack %s"
-                    % (stkResiduals, stkAutoCorrelations, fnAutoCorrelations), numberOfMpi=1)
-        self.runJob("xmipp_metadata_utilities", '-i %s --operate rename_column "image imageResidual"'
-                    % fnAutoCorrelations, numberOfMpi=1)
-        self.runJob("xmipp_metadata_utilities", '-i %s --set join %s imageResidual'
-                    % (anglesOutFn, fnAutoCorrelations), numberOfMpi=1)
-        cleanPath(fnAutoCorrelations)
-    
+        for volId, fnResiduals in self.fnResidualsDict.items():
+            # Evaluate each image
+            fnAutoCorrelations = self._getExtraPath("autocorrelations_%d.xmd" % volId)
+            stkAutoCorrelations = self._getExtraPath("autocorrelations_%d.mrcs" % volId)
+            stkResiduals = fnResiduals
+            anglesOutFn = self._getExtraPath("anglesCont_%d.xmd" % volId)
+            self.runJob("xmipp_image_residuals", " -i %s -o %s --save_metadata_stack %s"
+                        % (stkResiduals, stkAutoCorrelations, fnAutoCorrelations), numberOfMpi=1)
+            self.runJob("xmipp_metadata_utilities", '-i %s --operate rename_column "image imageResidual"'
+                        % fnAutoCorrelations, numberOfMpi=1)
+            self.runJob("xmipp_metadata_utilities", '-i %s --set join %s imageResidual'
+                        % (anglesOutFn, fnAutoCorrelations), numberOfMpi=1)
+            cleanPath(fnAutoCorrelations)
+
     def createOutputStep(self):
         fnImgs = self._getExtraPath('images.stk')
         if os.path.exists(fnImgs):
             cleanPath(fnImgs)
 
-        imgSet = self.inputSet.get()
-        imgFn = self._getExtraPath("anglesCont.xmd")
-        self.newAssignmentPerformed = os.path.exists(self._getExtraPath("angles.xmd"))
-        self.samplingRate = self.inputSet.get().getSamplingRate()
-
-        # Special case for 2D classes
-        if isinstance(imgSet, SetOfClasses2D):
-            outputSet = self._createSetOfClasses2D(imgSet.getImages())
-            outputSet.copyInfo(imgSet)
-            outputSet.appendFromClasses(imgSet, updateClassCallback=lambda clazz: self._updateClass(clazz, imgFn))
-
-        # Particles or Averages
-        else:
-            if isinstance(imgSet, SetOfAverages):
-                outputSet = self._createSetOfAverages()
+        imgSet = self.inputSet2D.get()
+        for volId, anglesFn in self.anglesDict.items():
+            imgFn = self._getExtraPath("anglesCont_%d.xmd" % volId)
+            self.newAssignmentPerformed = os.path.exists(anglesFn)
+            # Special case for 2D classes
+            if isinstance(imgSet, SetOfClasses2D):
+                outputSet = self._createSetOfClasses2D(imgSet.getImages(), suffix="_vol%d" % volId)
                 outputSet.copyInfo(imgSet)
+                outputSet.appendFromClasses(imgSet, updateClassCallback=lambda clazz: self._updateClass(clazz, imgFn))
+                self._classesInfo = dict()  # For every output you need to reset this value so the _updateClass works
+            # Particles or Averages
             else:
-                outputSet = self._createSetOfParticles()
+                if isinstance(imgSet, SetOfAverages):
+                    outputSet = self._createSetOfAverages(suffix="_vol%d" % volId)
+                else:
+                    outputSet = self._createSetOfParticles(suffix="_vol%d" % volId)
+                    if not self.newAssignmentPerformed:
+                        outputSet.setAlignmentProj()
+
                 outputSet.copyInfo(imgSet)
-                if not self.newAssignmentPerformed:
-                    outputSet.setAlignmentProj()
-            self.iterMd = md.iterRows(imgFn, md.MDL_ITEM_ID)
-            self.lastRow = next(self.iterMd)
-            outputSet.copyItems(imgSet, updateItemCallback=self._processRow)
+                outputSet.setDim((self.xDim, self.xDim, 1))
+                outputSet.setObjLabel("")  # To avoid renaming based on the label
+                self.iterMd = md.iterRows(imgFn, md.MDL_ITEM_ID)
+                self.lastRow = next(self.iterMd)
+                outputSet.copyItems(imgSet, updateItemCallback=self._processRow)
 
-        self._defineOutputs(reprojections=outputSet)
-        self._defineSourceRelation(self.inputSet, outputSet)
+            name = self.outputNamesDict[volId]
+            self._possibleOutputs[name] = outputSet
+            self._defineOutputs(**{name: outputSet})
+            self._defineSourceRelation(self.inputSet2D, outputSet)
+
+        if self.doRanking:
+            bestVolId = self.computeRankingVolumes(self._possibleOutputs)
+            if self.doExtraction:
+                volCl = self.inputSet3D.get().getItem("id", bestVolId)  # It may be a Volume or a 3D Class
+                outputParticles, outputVol = self._extractElementsFrom3D(volCl)
+                if outputParticles:
+                    particlesName = "particles_bestVol"
+                    self._defineOutputs(**{particlesName: outputParticles})
+                    self._defineSourceRelation(self.inputSet3D, outputParticles)
+
+                if outputVol:
+                    volName = "bestVolume"
+                    self._defineOutputs(**{volName: outputVol})
+                    self._store(outputVol)
+                    self._defineSourceRelation(self.inputSet3D, outputVol)
+
+        self.writeOutputDict()  # For visualization purpose
 
     def _updateClass(self, clazz, mdFile):
         """ Callback to update the class"""
@@ -253,7 +438,7 @@ def __setXmippImage(label):
 
     def _useProjMatching(self):
         """ Determine if it is necessary to perform projection matching step (because there is not input alignment)"""
-        imgSet = self.inputSet.get()
+        imgSet = self.inputSet2D.get()
         if not self.useAssignment or isinstance(imgSet, SetOfClasses2D) \
                 or (isinstance(imgSet, SetOfAverages) and not imgSet.hasAlignmentProj()) or \
                 (isinstance(imgSet, SetOfParticles) and not imgSet.hasAlignmentProj()):
@@ -277,29 +462,176 @@ def _computeResiduals(self, fnVol):
             cleanPath(fnNewParticles)
         if os.path.exists(mrcsresiduals):
             cleanPath(mrcsresiduals)
+        cleanPattern(self._getExtraPath("residual_part*.stk"))
+
+    def computeRankingVolumes(self, outputSetDict):
+        self.info('Ranking the best volumes...')
+        resultsDict = {}
+        meanCostDict = {}
+        for outname, outputSet in outputSetDict.items():
+            mdLabel = emlib.MDL_COST
+            xmippCostValues = {avg.getObjId(): getXmippAttribute(avg, mdLabel).get()
+                               for avg in outputSet}
+            resultsDict[outname] = xmippCostValues
+            meanCostDict[outname] = np.mean(list(xmippCostValues.values()))
+
+        # Find the key-value pair with the maximum value
+        bestVolume, bestScore = max(meanCostDict.items(), key=lambda item: item[1])
+        bestVolume = bestVolume.split('_')[-1]
+        volumeId = int(bestVolume[3:])
+        self.info(meanCostDict)
+        msg = "The volume with the best score has id %d and score %f" % (volumeId, bestScore)
+        self.info(msg)
+        self._storeSummaryInfo(msg)
+
+        return volumeId
+
+    def _extractElementsFrom3D(self, volCl):
+        """ Extract the elements (particles and/or volume) from the 3D class and create the output """
+        outputParticles, outputVol = self._getOutputSet()
+        if isinstance(volCl, Class3D):
+            self.info('The extraction 3D class have id %d with size %d' % (volCl.getObjId(), volCl.getSize()))
+            vol = volCl.getRepresentative().clone()
+            if outputParticles is not None:
+                # Go through all items and append them
+                for image in volCl:
+                    newImage = image.clone()
+                    outputParticles.append(newImage)
+        else:
+            self.info('The extraction 3D Volume have id %d' % volCl.getObjId())
+            vol = volCl.clone()
+
+        if outputVol is not None:
+            # Get the corresponding volume from the 3D class
+            outputVol.copyInfo(volCl)
+            outputVol.setLocation(vol.getLocation())
+            if vol.hasOrigin():
+                outputVol.setOrigin(vol.getOrigin())
+
+        return outputParticles, outputVol
+
+    def _getOutputSet(self):
+        """ Creates the output sets so they can be filled """
+        outputParticles = None
+        outputVol = None
+
+        if self.extractOption.get() == self.PARTICLES:
+            self.info("Creating set of particles")
+            outputParticles = createSetOfParticles(self.inputSet3D.get(), self._getPath())
+
+        elif self.extractOption.get() == self.VOLUME:
+            self.info("Creating volume")
+            outputVol = createRepresentativeVolume(self.inputSet3D.get())
+        else:  # Both
+            self.info("Creating both the volume and the set of particles")
+            outputParticles = createSetOfParticles(self.inputSet3D.get(), self._getPath())
+            outputVol = createRepresentativeVolume(self.inputSet3D.get())
+
+        return outputParticles, outputVol
+
+
+    def getOutputNamesDict(self):
+        return self.outputNamesDict
+
+    def writeOutputDict(self):
+        """Write a dictionary to a text file."""
+        dictionary = self.getOutputNamesDict()
+        filePath = self._getExtraPath(OUTPUTS_FN)
+        with open(filePath, 'w') as file:
+            json.dump(dictionary, file)
+
+    def readOutputDict(self):
+        """Read a dictionary from a text file."""
+        filePath = self._getExtraPath(OUTPUTS_FN)
+        with open(filePath, 'r') as file:
+            dictionary = json.load(file)
+        return dictionary
 
     # --------------------------- INFO functions --------------------------------------------
     def _summary(self):
         summary = []
-        summary.append("Images evaluated: %i" % self.inputSet.get().getSize())
-        summary.append("Volume: %s" % self.inputVolume.getNameId())
+        summary.append("Images evaluated: %i" % self.inputSet2D.get().getSize())
+        summary.append("Volume: %s" % self.inputSet3D.getNameId())
+
         return summary
-    
+
+    def _storeSummaryInfo(self, rankingMsg):
+        self.summaryVar.set(rankingMsg)
+
+    def _validate(self):
+        errors = []
+        if self.doRanking:
+            if not isinstance(self.inputSet3D.get(), EMSet):
+                errors.append("The input 3D must be a Set of Volumes or 3D classes to run the ranking option.")
+            if self.doExtraction:
+                if (isinstance(self.inputSet3D.get(), SetOfVolumes) and
+                        (self.extractOption.get() == self.PARTICLES or self.extractOption.get() == self.BOTH)):
+                    errors.append("The input 3D must be a Set of Classes 3D in order to extract its particles. "
+                                  "Please change to extract Volume option.")
+        return errors
+
     def _methods(self):
         methods = []
         if hasattr(self, 'outputParticles'):
             methods.append("We evaluated %i input images %s regarding to volume %s." %
-                           (self.inputSet.get().getSize(), self.getObjectTag('inputSet'),
-                            self.getObjectTag('inputVolume')))
+                           (self.inputSet2D.get().getSize(), self.getObjectTag('inputSet2D'),
+                            self.getObjectTag('inputSet3D')))
             methods.append("The residuals were evaluated according to their mean, variance and covariance structure "
                            "[Cherian2013].")
         return methods
     
     # --------------------------- UTILS functions --------------------------------------------
-    def _getDimensions(self):
-        imgSet = self.inputSet.get()
+    def _getDimensionsImages(self):
+        imgSet = self.inputSet2D.get()
         if isinstance(imgSet, SetOfClasses2D):
             xDim = imgSet.getImages().getDim()[0]
         else:
             xDim = imgSet.getDim()[0]
         return xDim
+
+    def _getDimensionsVol(self):
+        volSet = self.inputSet3D.get()
+        if isinstance(volSet, EMSet):
+            vol = volSet.getFirstItem()
+        else:
+            vol = volSet
+        xDimVol = vol.getDim()[0]
+        return xDimVol
+
+    def _getWorkingDimensions(self, srImages, srVol, doDownSample):
+        if doDownSample:
+            if srImages < srVol:
+                if NEW_SAMPLING_RATE < srVol:
+                    newDim = self._getDimensionsVol()
+                else:
+                    factor = srVol / NEW_SAMPLING_RATE
+                    newDim = self._getDimensionsVol() * factor
+            else:
+                if NEW_SAMPLING_RATE < srImages:
+                    newDim = self._getDimensionsImages()
+                else:
+                    factor = srImages / NEW_SAMPLING_RATE
+                    newDim = self._getDimensionsImages() * factor
+        else:
+            if srImages < srVol:
+                newDim = self._getDimensionsVol()
+            else:
+                newDim = self._getDimensionsImages()
+
+        return newDim
+
+#  ---------------------------- HELPERS --------------------------------------
+def createRepresentativeVolume(classesSet):
+    """ Creates a Volume from the corresponding set from the representative of a set of classes """
+    volInput = classesSet.getFirstItem()
+    vol = Volume()  # Create an instance of the volume
+    vol.copyInfo(volInput)
+    return vol
+
+
+def createSetOfParticles(classesSet, path):
+    """ Creates the corresponding set of particles from the input set of classes """
+    images = classesSet.getImages()
+    particles = SetOfParticles.create(outputPath=path)
+    particles.copyInfo(images)
+    return particles
\ No newline at end of file
diff --git a/xmipp3/protocols/protocol_convert_pdb.py b/xmipp3/protocols/protocol_convert_pdb.py
index 0697299ce..8fa0494da 100644
--- a/xmipp3/protocols/protocol_convert_pdb.py
+++ b/xmipp3/protocols/protocol_convert_pdb.py
@@ -40,12 +40,12 @@
 import pyworkflow.protocol.params as params
 import pyworkflow.protocol.constants as const
 from pyworkflow.utils import replaceBaseExt, removeExt, getExt, createLink, replaceExt, removeBaseExt
-from pyworkflow import UPDATED
+from pyworkflow import UPDATED, PROD
 
 class XmippProtConvertPdb(ProtInitialVolume):
     """ Convert atomic structure(s) into volume(s) """
     _label = 'convert pdbs to volumes'
-    _devStatus = UPDATED
+    _devStatus = PROD
     OUTPUT_NAME1 = "outputVolume"
     OUTPUT_NAME2 = "outputVolumes"
     OUTPUT_NAME3 = "outputPdb"
@@ -70,7 +70,7 @@ def _defineParams(self, form):
         This definition is also used to generate automatically the GUI.
         """
         # Defining condition string for x, y, z coords
-        coordsCondition = 'setSize and not vol'
+        coordsCondition = 'setSize and not vol and not isinstance(pdbObj, SetOfAtomStructs)'
 
         # Defining parallel arguments
         form.addParallelSection(threads=4, mpi=1)
@@ -254,6 +254,9 @@ def _validate(self):
         if self.numberOfMpi > 1:
             errors.append('MPI cannot be selected, because Scipion is going to drop support for it. Select threads instead.')
 
+        if isinstance(self.pdbObj.get(), SetOfAtomStructs) and not self.size.hasValue():
+            errors.append('Please set size when using SetOfAtomStructs as input')
+
         return errors
     
     # --------------------------- UTLIS functions --------------------------------------------
diff --git a/xmipp3/protocols/protocol_core_analysis.py b/xmipp3/protocols/protocol_core_analysis.py
index caeed6ba6..cd1a7f30a 100644
--- a/xmipp3/protocols/protocol_core_analysis.py
+++ b/xmipp3/protocols/protocol_core_analysis.py
@@ -30,8 +30,7 @@
 
 
 import pyworkflow.protocol.params as param
-import pyworkflow.protocol.constants as const
-from pyworkflow.utils.path import cleanPath, makePath
+from pyworkflow.utils.path import makePath
 
 import pwem.emlib.metadata as md
 from pwem.protocols import ProtClassify2D
@@ -91,11 +90,11 @@ def analyzeCore(self):
         inputMdName = join(fnLevel, 'level_classes.xmd')
         writeSetOfClasses2D(self.inputClasses.get(), inputMdName, writeParticles=True)
 
-        args = " --dir %s --root level --computeCore %f %f"%(self._getExtraPath(),
-                                                            self.thZscore, self.thPCAZscore)
+        args = " --dir %s --root level --computeCore %f %f" % (self._getExtraPath(),
+                                                               self.thZscore, self.thPCAZscore)
         self.runJob('xmipp_classify_CL2D_core_analysis', args)
         self.runJob("xmipp_classify_evaluate_classes", "-i %s"%\
-                    self._getExtraPath(join("level_00","level_classes_core.xmd")), numberOfMpi=1)
+                    self._getExtraPath(join("level_00", "level_classes_core.xmd")), numberOfMpi=1)
 
     #--------------------------- STEPS functions -------------------------------
     def _defineFileNames(self):
@@ -104,8 +103,8 @@ def _defineFileNames(self):
         myDict = {
                   'final_classes': self._getPath('classes2D%(sub)s.sqlite'),
                   'output_particles': self._getExtraPath('images.xmd'),
-                  'level_classes' : self.levelPath + 'level_classes%(sub)s.xmd',
-                  'level_images' : self.levelPath + 'level_images%(sub)s.xmd',
+                  'level_classes': self.levelPath + 'level_classes%(sub)s.xmd',
+                  'level_images': self.levelPath + 'level_images%(sub)s.xmd',
                   'classes_scipion': (self.levelPath + 'classes_scipion_level_'
                                                    '%(level)02d%(sub)s.sqlite'),
                   }
@@ -115,7 +114,7 @@ def createOutputStep(self):
         """ Store the SetOfClasses2D object
         resulting from the protocol execution.
         """
-        inputParticles = self.inputClasses.get().getImages()
+        inputParticles = self.inputClasses.get().getImagesPointer()
         level = self._getLastLevel()
         subset = CLASSES_CORE
 
@@ -143,7 +142,7 @@ def _methods(self):
         strline ='We calculated the class cores %s. [Sorzano2014]' % self.getObjectTag('outputClasses_core')
         return [strline]
     
-    #--------------------------- UTILS functions -------------------------------
+    # --------------------------- UTILS functions -------------------------------
     def _updateParticle(self, item, row):
         item.setClassId(row.getValue(md.MDL_REF))
         item.setTransform(rowToAlignment(row, ALIGN_2D))
diff --git a/xmipp3/protocols/protocol_ctf_micrographs.py b/xmipp3/protocols/protocol_ctf_micrographs.py
index fe1679f88..c83648414 100644
--- a/xmipp3/protocols/protocol_ctf_micrographs.py
+++ b/xmipp3/protocols/protocol_ctf_micrographs.py
@@ -132,7 +132,7 @@ def _defineProcessParams(self, form):
                            'and if it fails, -1.')
         form.addParam('refineAmplitudeContrast', params.BooleanParam, default=False,
                       label='Allow amplitude constrast refinement',
-                      help = 'The amplitude contrast is normally kept fixed, but in'
+                      help = 'The amplitude contrast is normally kept fixed, but in '
                              'some experiments it has been found that refining it '
                              'might result in some improvement in the final FSC. '
                              'This is not a standard practice, and should be used with caution')
@@ -289,14 +289,6 @@ def _getFn(key):
             sys.stderr.write("xmipp_ctf_estimate_from_micrograph has " \
                              "failed with micrograph %s" % finalName)
 
-    def _reEstimateCTF(self, mic, ctfModel):
-        """ Run the estimate CTF program """
-        self._prepareRecalCommand(ctfModel)
-        # CTF estimation with Xmipp
-        self.runJob(self._program, self._args % self._params)
-        mic = ctfModel.getMicrograph()
-        self.evaluateSingleMicrograph(mic)
-
     def _createOutputStep(self):
         pass
 
diff --git a/xmipp3/protocols/protocol_deep_center.py b/xmipp3/protocols/protocol_deep_center.py
new file mode 100644
index 000000000..16d72580e
--- /dev/null
+++ b/xmipp3/protocols/protocol_deep_center.py
@@ -0,0 +1,148 @@
+# **************************************************************************
+# *
+# * Authors:     COS Sorzano, Erney Ramirez and Adrian Sansinena
+# *
+# * 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'
+# *
+# **************************************************************************
+
+from pyworkflow import VERSION_3_0
+from pyworkflow.protocol.params import (PointerParam, StringParam, FloatParam,
+                                        IntParam, BooleanParam, GPU_LIST)
+from pyworkflow.protocol.constants import LEVEL_ADVANCED
+from pyworkflow.utils import Message
+from pyworkflow.utils.path import createLink
+from pwem.protocols import ProtAlign2D
+from pwem.objects import String
+from xmipp3.convert import readSetOfParticles, writeSetOfParticles
+import os
+from os import remove
+import xmipp3
+from xmipp_base import XmippScript
+from pyworkflow import BETA, UPDATED, NEW, PROD
+
+class XmippProtDeepCenter(ProtAlign2D, xmipp3.XmippProtocol):
+    """Center a set of particles in 2D using a neural network. The particles remain the same, but their alignment
+       includes an approximate shift to place them in the center."""
+    _lastUpdateVersion = VERSION_3_0
+    _conda_env = 'xmipp_DLTK_v1.0'
+    _label = 'deep center'
+    _devStatus = BETA
+
+    def __init__(self, **args):
+        ProtAlign2D.__init__(self, **args)
+
+    # --------------------------- DEFINE param functions --------------------------------------------
+    def _defineParams(self, form):
+        form.addParallelSection(threads=1, mpi=4)
+
+        form.addHidden(GPU_LIST, StringParam, default='0',
+                       expertLevel=LEVEL_ADVANCED,
+                       label="Choose GPU IDs",
+                       help="GPU may have several cores. Set it to zero"
+                            " if you do not know what we are talking about."
+                            " First core index is 0, second 1 and so on.")
+
+        form.addSection(label=Message.LABEL_INPUT)
+
+        form.addParam('inputParticles', PointerParam, label="Input images",
+                      pointerClass='SetOfParticles',
+                      help='The set does not need to be centered or have alignment parameters')
+
+        form.addParam('sigma', FloatParam,
+                      label="Shift sigma",
+                      default=5,
+                      help="In pixels. This is used to generate artificially shifted particles.")
+
+        form.addParam('precision', FloatParam,
+                      label="Precision",
+                      default=0.5,
+                      help="In pixels.")
+
+        form.addSection(label="Training")
+
+        form.addParam('trainSetSize', IntParam, label="Train set size", default=5000,
+                      help='How many particles to use for training. Set to -1 for all of them')
+
+        form.addParam('numEpochs', IntParam,
+                      label="Number of epochs",
+                      default=100,
+                      expertLevel=LEVEL_ADVANCED,
+                      help="Number of epochs for training.")
+
+        form.addParam('batchSize', IntParam,
+                      label="Batch size for training",
+                      default=32,
+                      expertLevel=LEVEL_ADVANCED,
+                      help="Batch size for training.")
+
+        form.addParam('learningRate', FloatParam,
+                      label="Learning rate",
+                      default=0.001,
+                      expertLevel=LEVEL_ADVANCED,
+                      help="Learning rate for training.")
+
+    # --------------------------- INSERT steps functions --------------------------------------------
+    def _insertAllSteps(self):
+        self.fnImgs = self._getTmpPath('imgs.xmd')
+        self.fnImgsTrain = self._getTmpPath('imgsTrain.xmd')
+        if self.useQueueForSteps() or self.useQueue():
+            myStr = os.environ["CUDA_VISIBLE_DEVICES"]
+        else:
+            myStr = self.gpuList.get()
+            os.environ["CUDA_VISIBLE_DEVICES"] = self.gpuList.get()
+        numGPU = myStr.split(',')
+        self._insertFunctionStep("convertInputStep", self.inputParticles.get())
+        self._insertFunctionStep("train", numGPU[0])
+        self._insertFunctionStep("predict", numGPU[0])
+        self._insertFunctionStep("createOutputStep")
+
+    # --------------------------- STEPS functions ---------------------------------------------------
+    def convertInputStep(self, inputSet):
+        writeSetOfParticles(inputSet, self.fnImgs)
+        if self.trainSetSize.get()>0:
+            self.runJob("xmipp_metadata_utilities","-i %s --operate random_subset %d -o %s"%\
+                        (self.fnImgs,self.trainSetSize, self.fnImgsTrain), numberOfMpi=1)
+        else:
+            createLink(self.fnImgs, self.fnImgsTrain)
+
+    def train(self, gpuId):
+        args = "-i %s --omodel %s --sigma %f --maxEpochs %d --batchSize %d --gpu %s --learningRate %f --precision %f"%\
+                (self.fnImgsTrain, self._getExtraPath("model.h5"), self.sigma, self.numEpochs, self.batchSize, gpuId,
+                 self.learningRate, self.precision)
+        self.runJob(f"xmipp_deep_center", args, numberOfMpi=1, env=self.getCondaEnv())
+
+    def predict(self, gpuId):
+        fnModel = self._getExtraPath("model.h5")
+        args = "-i %s --gpu %s --model %s -o %s" % (self.fnImgs, gpuId, fnModel, self._getExtraPath('particles.xmd'))
+        self.runJob("xmipp_deep_center_predict", args, numberOfMpi=1, env=self.getCondaEnv())
+
+    def createOutputStep(self):
+        fnPredict = self._getExtraPath("particles.xmd")
+        outputSet = self._createSetOfParticles()
+        readSetOfParticles(fnPredict, outputSet)
+        outputSet.copyInfo(self.inputParticles.get())
+        outputSet.setAlignment2D()
+        self._defineOutputs(outputParticles=outputSet)
+        self._store(outputSet)
+        self._defineSourceRelation(self.inputParticles.get(), outputSet)
+
+
diff --git a/xmipp3/protocols/protocol_deep_center_assignment.py b/xmipp3/protocols/protocol_deep_center_assignment.py
deleted file mode 100644
index 2cfe2d4bd..000000000
--- a/xmipp3/protocols/protocol_deep_center_assignment.py
+++ /dev/null
@@ -1,261 +0,0 @@
-# **************************************************************************
-# *
-# * Authors:     COS Sorzano, Erney Ramirez and Adrian Sansinena
-# *
-# * 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'
-# *
-# **************************************************************************
-
-from pyworkflow import VERSION_3_0
-from pyworkflow.protocol.params import (PointerParam, StringParam, FloatParam,
-                                        IntParam, BooleanParam, GPU_LIST)
-from pyworkflow.protocol.constants import LEVEL_ADVANCED
-from pyworkflow.utils import Message
-from pwem.protocols import ProtAlign2D
-from pwem.objects import String
-from xmipp3.convert import readSetOfParticles, writeSetOfParticles
-import os
-from os import remove
-import xmipp3
-from xmipp_base import XmippScript
-from pyworkflow import BETA, UPDATED, NEW, PROD
-
-class XmippProtDeepCenterAssignmentPredictBase(ProtAlign2D, xmipp3.XmippProtocol):
-    """Super protocol of Deep Center and Deep Global Assignment."""
-    _lastUpdateVersion = VERSION_3_0
-    _conda_env = 'xmipp_DLTK_v1.0'
-    _trainingResizedFileName = 'trainingResized'
-
-    def __init__(self, **args):
-        ProtAlign2D.__init__(self, **args)
-
-    # --------------------------- DEFINE param functions --------------------------------------------
-    def _defineParams(self, form):
-        form.addParallelSection(threads=1, mpi=1)
-
-        form.addHidden(GPU_LIST, StringParam, default='0',
-                       expertLevel=LEVEL_ADVANCED,
-                       label="Choose GPU IDs",
-                       help="GPU may have several cores. Set it to zero"
-                            " if you do not know what we are talking about."
-                            " First core index is 0, second 1 and so on.")
-
-        form.addSection(label=Message.LABEL_INPUT)
-
-        form.addParam('inputImageSet', PointerParam, label="Input Image set",
-                      pointerClass='SetOfParticles',
-                      help='The set of particles to predict')
-
-        form.addParam('trainModels', BooleanParam, label="Train models",
-                      pointerClass='SetOfParticles', default=False,
-                      help='Choose if you want to train a model using a centered set of particles')
-        
-        trainingGroup = form.addGroup('Training parameters', condition='trainModels==True')
-
-        trainingGroup.addParam('numModels', IntParam,
-                      label="Number of models", default=5,
-                      help="The maximum number of model available in xmipp is 5.")
-
-        trainingGroup.addParam('inputTrainSet', PointerParam, label="Input train set",
-                      pointerClass='SetOfParticles', allowsNull=True,
-                      pointerCondition='hasAlignment2D or hasAlignmentProj',
-                      help='The set of particles to train the models')
-
-        trainingGroup.addParam('numEpochs', IntParam,
-                      label="Number of epochs",
-                      default=25,
-                      expertLevel=LEVEL_ADVANCED,
-                      help="Number of epochs for training.")
-
-        trainingGroup.addParam('batchSize', IntParam,
-                      label="Batch size for training",
-                      default=32,
-                      expertLevel=LEVEL_ADVANCED,
-                      help="Batch size for training.")
-
-        trainingGroup.addParam('learningRate', FloatParam,
-                      label="Learning rate",
-                      default=0.001,
-                      expertLevel=LEVEL_ADVANCED,
-                      help="Learning rate for training.")
-
-        trainingGroup.addParam('sigma', FloatParam,
-                      label="Image shifting",
-                      default=5,
-                      expertLevel=LEVEL_ADVANCED,
-                      help="A measure of the number of pixels that particles can be shifted in each direction from the center.")
-
-        trainingGroup.addParam('patience', IntParam,
-                      label="Patience",
-                      default=5,
-                      expertLevel=LEVEL_ADVANCED,
-                      help="Training will be stopped if the number of epochs without improvement is greater than "
-                           "patience.")
-
-        form.addSection(label='Consensus')
-
-        form.addParam('tolerance', IntParam,
-                      label="Tolerance in pixels", default=3,
-                      help="Max difference between predictions and their mean value.")
-
-        form.addParam('maxModels', IntParam,
-                      label="Maximum number of models dropped per particle", default=0,
-                      help="If more models are dropped, the particle is discarded.")
-
-    # --------------------------- INSERT steps functions --------------------------------------------
-    def _insertAllSteps(self):
-        self.predictImgsFn = self._getExtraPath('predict_input_imgs.xmd')
-        if self.useQueueForSteps() or self.useQueue():
-            myStr = os.environ["CUDA_VISIBLE_DEVICES"]
-        else:
-            myStr = self.gpuList.get()
-            os.environ["CUDA_VISIBLE_DEVICES"] = self.gpuList.get()
-        numGPU = myStr.split(',')
-        self._insertFunctionStep("convertStep", self.inputImageSet.get())
-        self._insertFunctionStep("predict", numGPU[0], self.predictImgsFn)
-        self._insertFunctionStep('createOutputStep')
-    
-    def insertTrainSteps(self):
-        self.trainImgsFn = self._getExtraPath('train_input_imgs.xmd')
-        if self.useQueueForSteps() or self.useQueue():
-            myStr = os.environ["CUDA_VISIBLE_DEVICES"]
-        else:
-            myStr = self.gpuList.get()
-            os.environ["CUDA_VISIBLE_DEVICES"] = self.gpuList.get()
-        numGPU = myStr.split(',')
-
-        self._insertFunctionStep("convertTrainStep")
-        self._insertFunctionStep("train", numGPU[0])
-
-    # --------------------------- STEPS functions ---------------------------------------------------
-    def convertStep(self, inputSet):
-        self.Xdim = 128
-        writeSetOfParticles(inputSet, self.predictImgsFn)
-        self.runJob("xmipp_image_resize",
-                    "-i %s -o %s --save_metadata_stack %s --fourier %d" %
-                    (self.predictImgsFn,
-                     self._getExtraPath(f'{self._trainingResizedFileName}.stk'),
-                     self._getExtraPath(f'{self._trainingResizedFileName}.xmd'),
-                     self.Xdim), numberOfMpi=self.numberOfThreads.get() * self.numberOfMpi.get())
-    
-    def convertTrainStep(self):
-        self.convertStep(self.inputTrainSet.get())
-    
-    def train(self, gpuId, mode="", orderSymmetry=None):
-        args = "%s %s %f %d %d %s %d %f %d %d" % (
-            self._getExtraPath(f"{self._trainingResizedFileName}.xmd"), self._getExtraPath("model"), self.sigma.get(),
-            self.numEpochs, self.batchSize.get(), gpuId, self.numModels.get(), self.learningRate.get(),
-            self.patience.get(), 0)
-        if orderSymmetry:
-            args += " " + str(orderSymmetry)
-        self.runJob(f"xmipp_deep_{mode}", args, numberOfMpi=1, env=self.getCondaEnv())
-
-        remove(self._getExtraPath(f"{self._trainingResizedFileName}.xmd"))
-        remove(self._getExtraPath(f"{self._trainingResizedFileName}.stk"))
-
-    def predict(self, predictImgsFn, gpuId, mode="", inputModel="", trainedModel=True, orderSymmetry=None):
-        if mode != "center" or trainedModel:
-            args = "%s %s %s %s %s %d %d %d" % (
-                self._getExtraPath(f"{self._trainingResizedFileName}.xmd"), gpuId, self._getPath(), predictImgsFn,
-                inputModel, self.numModels.get(), self.tolerance.get(),
-                self.maxModels.get())
-            if orderSymmetry:
-                args += " " + str(orderSymmetry)
-        else:
-            args = "%s %s %s %s %s %d %d %d" % (
-                self._getExtraPath(f"{self._trainingResizedFileName}.xmd"), gpuId, self._getPath(), predictImgsFn,
-                os.path.join(XmippScript.getModel("deep_center"), 'modelCenter'), self.numModels.get(), self.tolerance.get(),
-                self.maxModels.get())
-        self.runJob(f"xmipp_deep_{mode}_predict", args, numberOfMpi=1, env=self.getCondaEnv())
-        remove(self._getExtraPath(f"{self._trainingResizedFileName}.xmd"))
-        remove(self._getExtraPath(f"{self._trainingResizedFileName}.stk"))
-
-    def createOutputStep(self):
-        imgFname = self._getPath('predict_results.xmd')
-        outputSet = self._createSetOfParticles()
-        readSetOfParticles(imgFname, outputSet)
-        outputSet.copyInfo(self.inputImageSet.get())
-        outputSet.setAlignmentProj()
-        self._defineOutputs(outputParticles=outputSet)
-        self._store(outputSet)
-        self._defineSourceRelation(self.inputImageSet.get(), outputSet)
-
-    # --------------------------- INFO functions --------------------------------
-    def _methods(self):
-        methods = []
-        if hasattr(self, 'outputParticles'):
-            methods.append("We learned a model to center particles from %i input images (%s)." \
-                           % (self.inputSet.get().getSize(), self.getObjectTag('inputSet')))
-        return methods
-
-class XmippProtDeepCenter(XmippProtDeepCenterAssignmentPredictBase):
-    """Predict the center particles using deep learning.""" 
-    _label = 'deep center'
-    _devStatus = BETA
-    # --------------------------- INSERT steps functions --------------------------------------------
-    def _insertAllSteps(self):
-        if self.trainModels.get():
-            self.insertTrainSteps()
-        super()._insertAllSteps()
-
-    # --------------------------- STEPS functions ---------------------------------------------------
-    def train(self, gpuId, mode="", orderSymmetry=None):
-        super().train(gpuId, mode="center")
-
-    def predict(self, predictImgsFn, gpuId, mode="", inputModel="", trainedModel=True, orderSymmetry=None):
-        super().predict(gpuId, predictImgsFn, mode="center",
-                                                           inputModel=self._getExtraPath("model"),
-                                                           trainedModel=self.trainModels.get())
-
-class XmippProtDeepGlobalAssignment(XmippProtDeepCenterAssignmentPredictBase):
-    """Predict Euler Angles using deep learning."""
-    _label = 'deep global assignment'
-    _devStatus = BETA
-
-    # --------------------------- DEFINE param functions --------------------------------------------
-    def _defineParams(self, form):
-        # Calling parent form
-        super()._defineParams(form)
-
-        # Adding extra conditions
-        trainModels = form.getParam('trainModels')
-        trainModels.condition = String('False')
-
-        # Always show training parameters group in this protocol
-        trainingGroup = form.getParam('Training_parameters')
-        trainingGroup.condition = String('True')
-
-        section = form.getSection(label=Message.LABEL_INPUT)
-        section.addParam('orderSymmetry', IntParam,
-                      label="Order of symmetry", default=1,
-                      help="Order of the group of the molecule.")
-    
-    # --------------------------- INSERT steps functions --------------------------------------------
-    def _insertAllSteps(self):
-        self.insertTrainSteps()
-        super()._insertAllSteps()
-
-    # --------------------------- STEPS functions ---------------------------------------------------
-    def train(self, gpuId, mode="", orderSymmetry=None):
-        super().train(gpuId, mode="global_assignment", orderSymmetry=self.orderSymmetry.get())
-
-    def predict(self, predictImgsFn, gpuId, mode="", inputModel="", trainedModel=True, orderSymmetry=None):
-        super().predict(gpuId, predictImgsFn, mode="global_assignment", inputModel=self._getExtraPath("model"), orderSymmetry=self.orderSymmetry.get())
diff --git a/xmipp3/protocols/protocol_denoise_particles.py b/xmipp3/protocols/protocol_denoise_particles.py
index 11cafb0a7..49c2938ec 100644
--- a/xmipp3/protocols/protocol_denoise_particles.py
+++ b/xmipp3/protocols/protocol_denoise_particles.py
@@ -80,25 +80,26 @@ def _insertAllSteps(self):
         """ Insert every step of the protocol"""
         
         # Convert input images if necessary
-        self._insertFunctionStep('denoiseImages', self.inputParticles.getObjId(), self.inputClasses.getObjId()) 
+        self._insertFunctionStep('denoiseImages')
         
         self._insertFunctionStep('createOutputStep')
 
     #--------------------------- STEPS functions --------------------------------------------
-    def denoiseImages(self, inputId, inputClassesId):
+    def denoiseImages(self):
         # We start preparing writing those elements we're using as input to keep them untouched
         imagesMd = self._getPath('images.xmd')
         writeSetOfParticles(self.inputParticles.get(), imagesMd)
         classesMd = self._getPath('classes.xmd')
         if isinstance(self.inputClasses.get(), SetOfAverages):
-            writeSetOfClasses2D(self.inputClasses.get(), classesMd)
+            writeSetOfParticles(self.inputClasses.get(), classesMd)
         else:
             writeSetOfClasses2D(self.inputClasses.get(), classesMd)
 
         fnRoot = self._getExtraPath('pca')
         fnRootDenoised = self._getExtraPath('imagesDenoised')
 
-        args = '-i Particles@%s --oroot %s --eigenvectors %d --maxImages %d' % (imagesMd, fnRoot, self.maxPCABases.get(), self.maxClasses.get())
+        args = ('-i Particles@%s --oroot %s --eigenvectors %d --maxImages %d'
+                % (imagesMd, fnRoot, self.maxPCABases.get(), self.maxClasses.get()))
         self.runJob("xmipp_image_rotational_pca", args)
 
         N=min(self.maxPCABases.get(), self.PCABases2Project.get())
diff --git a/xmipp3/protocols/protocol_extract_particles.py b/xmipp3/protocols/protocol_extract_particles.py
index 8fb8254c3..f226d9d41 100644
--- a/xmipp3/protocols/protocol_extract_particles.py
+++ b/xmipp3/protocols/protocol_extract_particles.py
@@ -322,7 +322,7 @@ def _validate(self):
             self.boxSize.set(self.getEven(self.boxSize))
 
         if self.doNormalize:
-            if self.backRadius > int(self.boxSize.get() / 2):
+            if self.backRadius >= int(self.boxSize.get() / 2):
                 errors.append("Background radius for normalization should be "
                               "equal or less than half of the box size.")
 
@@ -423,12 +423,13 @@ def _getNewSampling(self):
 
     def _getDownFactor(self):
         downFactor = 1
-        if self.resizeOption == self.RESIZE_FACTOR:
-            downFactor = self.downFactor.get()
-        elif self.resizeOption == self.RESIZE_SAMPLINGRATE:
-            downFactor = self.resizeSamplingRate.get()/self.getCoordSampling()
-        elif self.resizeOption == self.RESIZE_DIMENSIONS:
-            downFactor = self.boxSize.get()/self.resizeDim.get()
+        if self.doResize:
+            if self.resizeOption == self.RESIZE_FACTOR:
+                downFactor = self.downFactor.get()
+            elif self.resizeOption == self.RESIZE_SAMPLINGRATE:
+                downFactor = self.resizeSamplingRate.get()/self.getCoordSampling()
+            elif self.resizeOption == self.RESIZE_DIMENSIONS:
+                downFactor = self.boxSize.get()/self.resizeDim.get()
 
         return float(downFactor)
 
@@ -569,6 +570,7 @@ def readPartsFromMics(self, micList, outputParts):
                         # adding the variance and Gini coeff. value of the mic zone
                         setXmippAttributes(p, row, md.MDL_SCORE_BY_VAR)
                         setXmippAttributes(p, row, md.MDL_SCORE_BY_GINI)
+                        setXmippAttributes(p, row, md.MDL_LOCAL_AVERAGE)
                         if row.containsLabel(md.MDL_ZSCORE_DEEPLEARNING1):
                             setXmippAttributes(p, row, md.MDL_ZSCORE_DEEPLEARNING1)
 
diff --git a/xmipp3/protocols/protocol_generate_reprojections.py b/xmipp3/protocols/protocol_generate_reprojections.py
index b9b116127..82b56b150 100644
--- a/xmipp3/protocols/protocol_generate_reprojections.py
+++ b/xmipp3/protocols/protocol_generate_reprojections.py
@@ -28,7 +28,7 @@
 import os
 
 from pyworkflow import VERSION_2_0
-from pyworkflow.protocol.params import PointerParam
+from pyworkflow.protocol.params import PointerParam, BooleanParam
 from pyworkflow.utils.path import cleanPath
 from pwem.protocols import ProtAnalysis3D
 from pwem.objects import Image
@@ -60,6 +60,9 @@ def _defineParams(self, form):
         form.addParam('inputVolume', PointerParam, label="Volume to compare images to", important=True,
                       pointerClass='Volume',
                       help='Volume to be used for class comparison')
+        form.addParam('ignoreCTF', BooleanParam, default=True, label='Ignore CTF',
+                      help='By ignoring the CTF you will create projections more similar to what a person expects, '
+                           'while by using the CTF you will create projections more similar to what the microscope sees')
         form.addParallelSection(threads=0, mpi=8)
     
     #--------------------------- INSERT steps functions --------------------------------------------
@@ -68,16 +71,16 @@ def _insertAllSteps(self):
         self.imgsFn = self._getExtraPath('input_imgs.xmd')
         vol = self.inputVolume.get()
         
-        self._insertFunctionStep("convertStep", self.imgsFn)
+        self._insertFunctionStep("convertStep")
         imgSet = self.inputSet.get()
-        anglesFn=self.imgsFn
-        self._insertFunctionStep("produceProjections", vol.getFileName(),
+        anglesFn = self.imgsFn
+        self._insertFunctionStep("produceProjections",
                                  anglesFn,
                                  vol.getSamplingRate())
         self._insertFunctionStep("createOutputStep")
 
     #--------------------------- STEPS functions ---------------------------------------------------
-    def convertStep(self, imgsFn):
+    def convertStep(self):
         from xmipp3.convert import writeSetOfParticles
         imgSet = self.inputSet.get()
         writeSetOfParticles(imgSet, self.imgsFn)
@@ -86,21 +89,22 @@ def convertStep(self, imgsFn):
         img = ImageHandler()
         fnVol = self._getTmpPath("volume.vol")
         img.convert(self.inputVolume.get(), fnVol)
-        xdim=self.inputVolume.get().getDim()[0]
+        xdim = self.inputVolume.get().getDim()[0]
 
-        imgXdim=imgSet.getDim()[0]
-        if xdim!=imgXdim:
-            self.runJob("xmipp_image_resize","-i %s --dim %d"%(fnVol,imgXdim),numberOfMpi=1)
+        imgXdim = imgSet.getDim()[0]
+        if xdim != imgXdim:
+            self.runJob("xmipp_image_resize", "-i %s --dim %d" % (fnVol, imgXdim), numberOfMpi=1)
     
-    def produceProjections(self, fnVol, fnAngles, Ts):
+    def produceProjections(self, fnAngles, Ts):
         fnVol = self._getTmpPath("volume.vol")
-        anglesOutFn=self._getExtraPath("anglesCont.stk")
-        projectionsOutFn=self._getExtraPath("projections.stk")
-        args="-i %s -o %s --ref %s --oprojections %s --sampling %f " \
-             "--max_angular_change 90"%(fnAngles,anglesOutFn,fnVol,
-                                  projectionsOutFn,Ts)
+        anglesOutFn = self._getExtraPath("anglesCont.stk")
+        projectionsOutFn = self._getExtraPath("projections.stk")
+        args = "-i %s -o %s --ref %s --oprojections %s --sampling %f " %\
+             (fnAngles, anglesOutFn, fnVol, projectionsOutFn, Ts)
+        if self.ignoreCTF:
+            args += " --ignoreCTF "
         self.runJob("xmipp_angular_continuous_assign2", args)
-        fnNewParticles=self._getExtraPath("images.stk")
+        fnNewParticles = self._getExtraPath("images.stk")
         if os.path.exists(fnNewParticles):
             cleanPath(fnNewParticles)
     
diff --git a/xmipp3/protocols/protocol_movie_alignment_consensus.py b/xmipp3/protocols/protocol_movie_alignment_consensus.py
new file mode 100644
index 000000000..07712c924
--- /dev/null
+++ b/xmipp3/protocols/protocol_movie_alignment_consensus.py
@@ -0,0 +1,602 @@
+# **************************************************************************
+# *
+# * Authors:    Carlos Oscar Sorzano (coss@cnb.csic.es)
+# *             Daniel Marchán Torres (da.marchan@cnb.csic.es)  -- streaming version
+# *
+# * 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 'coss@cnb.csic.es'
+# *
+# **************************************************************************
+"""
+Consensus alignment protocol
+"""
+
+import os
+from datetime import datetime
+from pyworkflow.gui.plotter import Plotter
+import numpy as np
+from math import ceil
+try:
+    from itertools import izip
+except ImportError:
+    izip = zip
+from pwem.objects import SetOfMovies, SetOfMicrographs, MovieAlignment, Image
+from pyworkflow.object import Set
+import pyworkflow.protocol.params as params
+from pyworkflow.protocol import STEPS_PARALLEL, Protocol
+import pyworkflow.utils as pwutils
+from pwem.protocols import ProtAlignMovies
+from pyworkflow.protocol.constants import (STATUS_NEW)
+from xmipp3.convert import getScipionObj
+from pwem.constants import ALIGN_NONE
+from pyworkflow import BETA, UPDATED, NEW, PROD
+
+ACCEPTED = 'Accepted'
+DISCARDED = 'Discarded'
+
+
+class XmippProtConsensusMovieAlignment(ProtAlignMovies, Protocol):
+    """
+    Protocol to estimate the agreement between different movie alignment
+    algorithms in the Global Shifts.
+    """
+
+    _label = 'movie alignment consensus'
+    outputName = 'consensusAlignments'
+    _devStatus = NEW
+
+    def __init__(self, **args):
+        ProtAlignMovies.__init__(self, **args)
+        self.stepsExecutionMode = STEPS_PARALLEL
+
+    def _defineParams(self, form):
+        form.addSection(label='Input Consensus')
+        form.addParam('inputMovies1', params.PointerParam, pointerClass='SetOfMovies',
+                      label="Reference Aligned Movies", important=True,
+                      help='Select the aligned movies to evaluate (this first set will give the global shifts)')
+
+        form.addParam('inputMovies2', params.PointerParam,
+                      pointerClass='SetOfMovies',
+                      label="Secondary Aligned Movies",
+                      help='Shift to be compared with reference alignment')
+
+        form.addParam('minConsCorrelation', params.FloatParam, default=-1,
+                      label='Minimum consensus shifts correlation',
+                      help="Minimum value for the consensus correlations between shifts trajectories."
+                           "\nIf there are noticeable discrepancies between the two estimations below this correlation,"
+                           " it will be discarded. If this value is set to -1 no movies will be discarded."
+                           "\n Reasonable values are from 0.5 to 1 meaning consistency between estimations.")
+
+        form.addParam('trajectoryPlot', params.BooleanParam, default=False,
+                      label='Global Alignment Trajectory Plot',
+                      help="This will generate a plot for each movie where the reference and the secondary trajectory"
+                           "will be plot in the same graph with its correlation value.")
+
+        form.addParallelSection(threads=4, mpi=1)
+
+# --------------------------- INSERT steps functions -------------------------
+    def _insertAllSteps(self):
+        self.initializeParams()
+        movieSteps = self._insertNewMovieSteps(self.allMovies1.keys(),
+                                               self.allMovies2.keys(),
+                                               self.insertedDict)
+        self._insertFunctionStep('createOutputStep',
+                                 prerequisites=movieSteps, wait=True)
+
+    def createOutputStep(self):
+        self._closeOutputSet()
+
+    def initializeParams(self):
+        self.finished = False
+        self.insertedDict = {}
+        self.processedDict = []
+        self.movieFn1 = self.inputMovies1.get().getFileName()
+        self.movieFn2 = self.inputMovies2.get().getFileName()
+        self.micsFn = self._getMicsPath()
+        self.stats = {}
+        self.isStreamClosed = self.inputMovies1.get().isStreamClosed() and \
+                              self.inputMovies2.get().isStreamClosed()
+        self.samplingRate = self.inputMovies1.get().getSamplingRate()
+        self.acquisition = self.inputMovies1.get().getAcquisition()
+        self.allMovies1 = {movie.getObjId(): movie.clone() for movie
+                           in self._loadInputMovieSet(self.movieFn1).iterItems()}
+        self.allMovies2 = {movie.getObjId(): movie.clone() for movie
+                           in self._loadInputMovieSet(self.movieFn2).iterItems()}
+        pwutils.makePath(self._getExtraPath('DONE'))
+
+    def _getFirstJoinStepName(self):
+        # This function will be used for streaming, to check which is
+        # the first function that need to wait for all movies
+        # to have completed, this can be overriden in subclasses
+        # (e.g., in Xmipp 'sortPSDStep')
+        return 'createOutputStep'
+
+    def _getFirstJoinStep(self):
+        for s in self._steps:
+            if s.funcName == self._getFirstJoinStepName():
+                return s
+        return None
+
+    def _stepsCheck(self):
+        self._checkNewInput()
+        self._checkNewOutput()
+
+    def _checkNewInput(self):
+        # Check if there are new movies to process from the input set
+        self.lastCheck = getattr(self, 'lastCheck', datetime.now())
+        mTime = max(datetime.fromtimestamp(os.path.getmtime(self.movieFn1)),
+                    datetime.fromtimestamp(os.path.getmtime(self.movieFn2)))
+
+        self.debug('Last check: %s, modification: %s'
+                   % (pwutils.prettyTime(self.lastCheck),
+                      pwutils.prettyTime(mTime)))
+        # If the input movies.sqlite have not changed since our last check,
+        # it does not make sense to check for new input data
+        if self.lastCheck is not None and self.lastCheck > mTime:
+            return None
+
+        movieSet1 = self._loadInputMovieSet(self.movieFn1)
+        movieSet2 = self._loadInputMovieSet(self.movieFn2)
+
+        movieDict1 = {movie.getObjId(): movie.clone() for movie in movieSet1.iterItems()}
+        movieDict2 = {movie.getObjId(): movie.clone() for movie in movieSet2.iterItems()}
+
+        newIds1 = [idMovie for idMovie in movieDict1.keys() if idMovie not in self.processedDict]
+        self.allMovies1.update(movieDict1)
+
+        newIds2 = [idMovie for idMovie in movieDict2.keys() if idMovie not in self.processedDict]
+        self.allMovies2.update(movieDict2)
+
+        self.lastCheck = datetime.now()
+        self.isStreamClosed = movieSet1.isStreamClosed() and \
+                              movieSet2.isStreamClosed()
+
+        movieSet1.close()
+        movieSet2.close()
+
+        outputStep = self._getFirstJoinStep()
+
+        if len(set(self.allMovies1)) > len(set(self.processedDict)) and \
+           len(set(self.allMovies2)) > len(set(self.processedDict)):
+
+            fDeps = self._insertNewMovieSteps(newIds1, newIds2, self.insertedDict)
+
+            if outputStep is not None:
+                outputStep.addPrerequisites(*fDeps)
+
+            self.updateSteps()
+
+
+    def _insertNewMovieSteps(self, movies1Dict, movies2Dict, insDict):
+        deps = []
+
+        newIDs = list(set(movies1Dict).intersection(set(movies2Dict)))
+
+        for movieID in newIDs:
+            if movieID not in insDict:
+                stepId = self._insertFunctionStep('alignmentCorrelationMovieStep', movieID, prerequisites=[])
+                deps.append(stepId)
+                insDict[movieID] = stepId
+                self.processedDict.append(movieID)
+
+        return deps
+
+    def alignmentCorrelationMovieStep(self, movieId):
+        movie1 = self.allMovies1.get(movieId)
+        movie2 = self.allMovies2.get(movieId)
+        doneFn = self._getMovieDone(movieId)
+
+        if self.isContinued() and self._isMovieDone(movieId):
+            self.info("Skipping movie with ID: %s, seems to be done" % movieId)
+            return
+
+        # Clean old finished files
+        pwutils.cleanPath(doneFn)
+
+        if (movie1 is None) or (movie2 is None):
+            self.info('AlignmentCorrelationMovieStep movie1 or movie2 are None')
+            return
+
+        alignment1 = movie1.getAlignment()
+        alignment2 = movie2.getAlignment()
+        shiftX_1, shiftY_1 = alignment1.getShifts()
+        shiftX_2, shiftY_2 = alignment2.getShifts()
+
+        # Transformation of the shifts to calculate the shifts trajectory correlation
+        S1 = np.ones([3, len(shiftX_1)])
+        S2 = np.ones([3, len(shiftX_2)])
+
+        S1[0, :] = shiftX_1
+        S1[1, :] = shiftY_1
+        S2[0, :] = shiftX_2
+        S2[2, :] = shiftY_2
+
+        A = np.dot(np.dot(S1, S2.T), np.linalg.inv(np.dot(S2, S2.T)))
+        S2_p = np.dot(A, S2)
+
+        S1_cart = np.array([S1[0, :]/S1[2, :], S1[1, :]/S1[2, :]])
+        S2_p_cart = np.array([S2_p[0, :] / S2_p[2, :], S2_p[1, :] / S2_p[2, :]])
+        rmse_cart = np.sqrt((np.square(S1_cart - S2_p_cart)).mean())
+        maxe_cart = np.max(S1_cart - S2_p_cart)
+        corrX_cart = np.corrcoef(S1_cart[0, :], S2_p_cart[0, :])[0, 1]
+        corrY_cart = np.corrcoef(S1_cart[1, :], S2_p_cart[1, :])[0, 1]
+        corr_cart = np.min([corrY_cart, corrX_cart])
+
+        self.info('Root Mean Squared Error %f' % rmse_cart)
+        self.info('General Corr min(corrX, corrY) %f' % corr_cart)
+
+        if corr_cart >= self.minConsCorrelation.get():
+            self.info('Movie with id %d has a correlated alignment shift trajectory' %movieId)
+            fn = self._getMovieSelecFileAccepted()
+            with open(fn, 'a') as f:
+                f.write('%d T\n' % movieId)
+
+        elif corr_cart < self.minConsCorrelation.get():
+            self.info('Movie with id %d has discrepancy in the alignment with correlation %f' % (movieId, corr_cart))
+            fn = self._getMovieSelecFileDiscarded()
+            with open(fn, 'a') as f:
+                f.write('%d F\n' % movieId)
+
+        stats_loc = {'shift_corr': corr_cart, 'shift_corr_X': corrX_cart, 'shift_corr_Y': corrY_cart,
+                     'max_error': maxe_cart, 'rmse_error': rmse_cart, 'S1_cart': S1_cart, 'S2_p_cart': S2_p_cart}
+
+        self.stats[movieId] = stats_loc
+        self._store()
+        # Mark this movie as finished
+        open(doneFn, 'w').close()
+
+    def _checkNewOutput(self):
+        """ Check for already selected movies and update the output set. """
+        # Load previously done items (from text file)
+        doneListDiscarded = self._readCertainDoneList(DISCARDED)
+        doneListAccepted = self._readCertainDoneList(ACCEPTED)
+        # Check for newly done items
+        movieListIdAccepted = self._readtMovieId(True)
+        movieListIdDiscarded = self._readtMovieId(False)
+
+        newDoneAccepted = [movieId for movieId in movieListIdAccepted
+                           if movieId not in doneListAccepted]
+        newDoneDiscarded = [movieId for movieId in movieListIdDiscarded
+                            if movieId not in doneListDiscarded]
+
+        firstTimeAccepted = len(doneListAccepted) == 0
+        firstTimeDiscarded = len(doneListDiscarded) == 0
+
+        allDone = len(doneListAccepted) + len(doneListDiscarded) +\
+                  len(newDoneAccepted) + len(newDoneDiscarded)
+
+        # We have finished when there is not more input movies (stream closed)
+        # and the number of processed movies is equal to the number of inputs
+        maxMovieSize = min(len(self.allMovies1), len(self.allMovies2))
+        self.finished = (self.isStreamClosed and allDone == maxMovieSize)
+        streamMode = Set.STREAM_CLOSED if self.finished else Set.STREAM_OPEN
+
+        def readOrCreateOutputs(doneList, newDone, label=''):
+            if len(doneList) > 0 or len(newDone) > 0:
+                with self._lock:
+                    movSet = self._loadOutputSet(SetOfMovies, 'movies'+label+'.sqlite')
+                    micSet = self._loadOutputSet(SetOfMicrographs, 'micrographs'+label+'.sqlite')
+                    label = ACCEPTED if label == '' else DISCARDED
+                    self.fillOutput(movSet, micSet, newDone, label)
+                    movSet.setSamplingRate(self.samplingRate)
+                    micSet.setSamplingRate(self.samplingRate)
+                    micSet.setAcquisition(self.acquisition.clone())
+                    movSet.setAcquisition(self.acquisition.clone())
+
+                return movSet, micSet
+            return None, None
+
+        movieSet, micSet = readOrCreateOutputs(doneListAccepted, newDoneAccepted)
+        movieSetDiscarded, micSetDiscarded = readOrCreateOutputs(doneListDiscarded, newDoneDiscarded, DISCARDED)
+
+        if not self.finished and not newDoneDiscarded and not newDoneAccepted:
+        # If we are not finished and no new output have been produced
+        # it does not make sense to proceed and updated the outputs
+        # so we exit from the function here
+            return
+
+        def updateRelationsAndClose(movieSet, micSet, first, label=''):
+            if os.path.exists(self._getPath('movies'+label+'.sqlite')):
+                micsAttrName = 'outputMicrographs'+label
+                self._updateOutputSet(micsAttrName, micSet, streamMode)
+                self._updateOutputSet('outputMovies'+label, movieSet, streamMode)
+
+                if first:
+                    # We consider that Movies are 'transformed' into the Micrographs
+                    # This will allow to extend the micrograph associated to a set of
+                    # movies to another set of micrographs generated from a
+                    # different movie alignment
+                    self._defineTransformRelation(self.inputMovies1, micSet)
+
+                micSet.close()
+                movieSet.close()
+
+        updateRelationsAndClose(movieSet, micSet, firstTimeAccepted)
+        updateRelationsAndClose(movieSetDiscarded, micSetDiscarded, firstTimeDiscarded, DISCARDED)
+
+        if self.finished:  # Unlock createOutputStep if finished all jobs
+            outputStep = self._getFirstJoinStep()
+            if outputStep and outputStep.isWaiting():
+                outputStep.setStatus(STATUS_NEW)
+
+    def fillOutput(self, movieSet, micSet, newDone, label):
+        if newDone:
+            inputMovieSet = self._loadInputMovieSet(self.movieFn1)
+            inputMicSet = self._loadInputMicrographSet(self.micsFn)
+
+            for movieId in newDone:
+                movie = inputMovieSet[movieId].clone()
+                mic = inputMicSet[movieId].clone()
+
+                movie.setEnabled(self._getEnable(movieId))
+                mic.setEnabled(self._getEnable(movieId))
+                alignment1 = movie.getAlignment()
+                shiftX_1, shiftY_1 = alignment1.getShifts()
+                setAttribute(mic, '_alignment_corr', self.stats[movieId]['shift_corr'])
+                setAttribute(mic, '_alignment_rmse_error', self.stats[movieId]['rmse_error'])
+                setAttribute(mic, '_alignment_max_error', self.stats[movieId]['max_error'])
+                alignment = MovieAlignment(xshifts=shiftX_1, yshifts=shiftY_1)
+                movie.setAlignment(alignment)
+
+                self._writeCertainDoneList(movieId, label)
+
+                if self.trajectoryPlot.get():
+                    firstFrame, _, _ = self.inputMovies1.get().getFramesRange()
+                    self._createAndSaveTrajectoriesPlot(movieId, firstFrame, self.samplingRate)
+                    mic.plotCart = Image()
+                    mic.plotCart.setFileName(self._getTrajectoriesPlot(movieId))
+
+                movieSet.append(movie)
+                micSet.append(mic)
+
+            inputMovieSet.close()
+            inputMicSet.close()
+
+    def _loadOutputSet(self, SetClass, baseName, fixSampling=True):
+        """
+        Load the output set if it exists or create a new one.
+        """
+        setFile = self._getPath(baseName)
+
+        if os.path.exists(setFile) and os.path.getsize(setFile) > 0:
+            outputSet = SetClass(filename=setFile)
+            outputSet.loadAllProperties()
+            outputSet.enableAppend()
+        else:
+            outputSet = SetClass(filename=setFile)
+            outputSet.setStreamState(outputSet.STREAM_OPEN)
+
+            inputMovies = self.inputMovies1.get()
+            outputSet.copyInfo(inputMovies)
+
+            if fixSampling:
+                newSampling = inputMovies.getSamplingRate() * self._getBinFactor()
+                outputSet.setSamplingRate(newSampling)
+
+        return outputSet
+
+    def _loadInputMovieSet(self, moviesFn):
+        self.debug("Loading input db: %s" % moviesFn)
+        movieSet = SetOfMovies(filename=moviesFn)
+        movieSet.loadAllProperties()
+        movieSet.close()
+        self.debug("Closed db.")
+        return movieSet
+
+    def _loadInputMicrographSet(self, micsFn):
+        self.debug("Loading input db: %s" % micsFn)
+        micSet = SetOfMicrographs(filename=micsFn)
+        micSet.loadAllProperties()
+        micSet.close()
+        self.debug("Closed db.")
+        return micSet
+
+    def _summary(self):
+        # return message
+        pass
+
+    def _validate(self):
+        """ The function of this hook is to add some validation before the
+        protocol is launched to be executed. It should return a list of
+        errors. If the list is empty the protocol can be executed.
+        """
+        errors = []
+        if (self.inputMovies1.get().hasAlignment() == ALIGN_NONE) or \
+           (self.inputMovies2.get().hasAlignment() == ALIGN_NONE):
+            errors.append("The inputs ( _Input Movies 1_ or _Input Movies 2_ must be aligned before")
+
+        return errors
+
+
+    # ------------------------------------ Utils functions ------------------------------------
+    def _isMovieDone(self, id):
+        """ A movie is done if the marker file exists. """
+        return os.path.exists(self._getMovieDone(id))
+
+    def _getMovieDone(self, id):
+        """ Return the file that is used as a flag of termination. """
+        return self._getExtraPath('DONE', 'movie_%06d.TXT' % id)
+
+    def _readDoneList(self):
+        """ Read from a file the id's of the items that have been done. """
+        doneFile = self._getAllDone()
+        doneList = []
+        # Check what items have been previously done
+        if os.path.exists(doneFile):
+            with open(doneFile) as f:
+                doneList += [int(line.strip()) for line in f]
+        return doneList
+
+    def _getAllDone(self):
+        return self._getExtraPath('DONE_all.TXT')
+
+    def _writeDoneList(self, partList):
+        """ Write to a text file the items that have been done. """
+        with open(self._getAllDone(), 'a') as f:
+            for part in partList:
+                f.write('%d\n' % part.getObjId())
+
+    def _getMicsPath(self):
+        prot1 = self.inputMovies1.getObjValue()  # pointer to previous protocol
+
+        if hasattr(prot1, 'outputMicrographs'):
+            path1 = prot1.outputMicrographs.getFileName()
+            if os.path.getsize(path1) > 0:
+                return path1
+        elif hasattr(prot1, 'outputMicrographsDoseWeighted'):
+            path2 = prot1.outputMicrographsDoseWeighted.getFileName()
+            if os.path.getsize(path2) > 0:
+                return path2
+        else:
+            return None
+
+    def _readCertainDoneList(self, label):
+        """ Read from a text file the id's of the items
+        that have been done. """
+        doneFile = self._getCertainDone(label)
+        doneList = []
+        # Check what items have been previously done
+        if os.path.exists(doneFile):
+            with open(doneFile) as f:
+                doneList += [int(line.strip()) for line in f]
+        return doneList
+
+    def _writeCertainDoneList(self, movieId, label):
+        """ Write to a text file the items that have been done. """
+        doneFile = self._getCertainDone(label)
+        with open(doneFile, 'a') as f:
+            f.write('%d\n' % movieId)
+
+    def _createAndSaveTrajectoriesPlot(self, movieId, first, pixSize):
+        """ Write to a text file the items that have been done. """
+        stats = self.stats[movieId]
+        fn = self._getExtraPath('global_trajectories_%d' %movieId+'_plot_cart.png')
+        shift_X1 = stats['S1_cart'][0, :]
+        shift_Y1 = stats['S1_cart'][1, :]
+        shift_X2 = stats['S2_p_cart'][0, :]
+        shift_Y2 = stats['S2_p_cart'][1, :]
+        # ---------------- PLOT -----------------------
+        sumMeanX1 = []
+        sumMeanY1= []
+        sumMeanX2 = []
+        sumMeanY2 = []
+
+        def px_to_ang(px):
+            y1, y2 = px.get_ylim()
+            x1, x2 = px.get_xlim()
+            ax_ang2.set_ylim(y1 * pixSize, y2 * pixSize)
+            ax_ang.set_xlim(x1 * pixSize, x2 * pixSize)
+            ax_ang.figure.canvas.draw()
+            ax_ang2.figure.canvas.draw()
+
+        figureSize = (6, 4)
+        plotter = Plotter(*figureSize)
+        figure = plotter.getFigure()
+        ax_px = figure.add_subplot(111)
+        ax_px.grid()
+
+        ax_px.set_xlabel('Shift x (px)')
+        ax_px.set_ylabel('Shift y (px)')
+
+        ax_px.set_xlabel('Shift x (px) (CorrX:%.3f)' % stats['shift_corr_X'])
+        ax_px.set_ylabel('Shift y (px) (CorrX:%.3f)' % stats['shift_corr_Y'])
+
+        ax_ang = ax_px.twiny()
+        ax_ang.set_xlabel('Shift x (A)')
+        ax_ang2 = ax_px.twinx()
+        ax_ang2.set_ylabel('Shift y (A)')
+
+        i = first
+        # The output and log files list the shifts relative to the first frame.
+        # ROB unit seems to be pixels since sampling rate is only asked
+        # by the program if dose filtering is required
+        skipLabels = ceil(len(shift_X1) / 10.0)
+        labelTick = 1
+
+        for x1, y1, x2, y2 in zip(shift_X1, shift_Y1, shift_X2, shift_Y2):
+            sumMeanX1.append(x1)
+            sumMeanY1.append(y1)
+            sumMeanX2.append(x2)
+            sumMeanY2.append(y2)
+
+            if labelTick == 1:
+                ax_px.text(x1 - 0.02, y1 + 0.02, str(i))
+                labelTick = skipLabels
+            else:
+                labelTick -= 1
+            i += 1
+
+        # automatically update lim of ax_ang when lim of ax_px changes.
+        ax_px.callbacks.connect("ylim_changed", px_to_ang)
+        ax_px.callbacks.connect("xlim_changed", px_to_ang)
+
+        ax_px.plot(sumMeanX1, sumMeanY1, color='b', label='reference shifts')
+        ax_px.plot(sumMeanX2, sumMeanY2, color='r', label='target shifts')
+        ax_px.plot(sumMeanX1, sumMeanY1, 'yo')
+        ax_px.plot(sumMeanX1[0], sumMeanY1[0], 'ro', markersize=10, linewidth=0.5)
+        ax_px.set_title('Global frame alignment')
+
+        ax_px.legend()
+        plotter.tightLayout()
+        plotter.savefig(fn)
+        plotter.close()
+
+    def _getTrajectoriesPlot(self, movieId):
+        """ Write to a text file the items that have been done. """
+        return self._getExtraPath('global_trajectories_%d' %movieId+'_plot_cart.png')
+
+    def _getCertainDone(self, label):
+        return self._getExtraPath('DONE_'+label+'.TXT')
+
+    def _getMovieSelecFileAccepted(self):
+        return self._getExtraPath('selection-movie-accepted.txt')
+
+    def _getMovieSelecFileDiscarded(self):
+        return self._getExtraPath('selection-movie-discarded.txt')
+
+    def _readtMovieId(self, accepted):
+        if accepted:
+            fn = self._getMovieSelecFileAccepted()
+        else:
+            fn = self._getMovieSelecFileDiscarded()
+        moviesList = []
+        # Check what items have been previously done
+        if os.path.exists(fn):
+            with open(fn) as f:
+                moviesList += [int(line.strip().split()[0]) for line in f]
+        return moviesList
+
+    def _getEnable(self, movieId):
+        fn = self._getMovieSelecFileAccepted()
+        # Check what items have been previously done
+        if os.path.exists(fn):
+            with open(fn) as f:
+                for line in f:
+                    if movieId == int(line.strip().split()[0]):
+                        if line.strip().split()[1] == 'T':
+                            return True
+                        else:
+                            return False
+
+def setAttribute(obj, label, value):
+    if value is None:
+        return
+    setattr(obj, label, getScipionObj(value))
+
diff --git a/xmipp3/protocols/protocol_movie_dose_analysis.py b/xmipp3/protocols/protocol_movie_dose_analysis.py
index 200190c37..60c01e65e 100644
--- a/xmipp3/protocols/protocol_movie_dose_analysis.py
+++ b/xmipp3/protocols/protocol_movie_dose_analysis.py
@@ -27,9 +27,9 @@
 import matplotlib.pyplot as plt
 import numpy as np
 import os
-from pyworkflow import VERSION_3_0, NEW
+from pyworkflow import VERSION_3_0
 from pyworkflow.object import Set
-from pyworkflow.protocol import STEPS_PARALLEL
+from pyworkflow.protocol import Protocol
 from pyworkflow.protocol.params import (PointerParam, IntParam, FloatParam, LEVEL_ADVANCED)
 from pyworkflow.utils.properties import Message
 import pyworkflow.protocol.constants as cons
@@ -37,21 +37,18 @@
 from pwem.objects import SetOfMovies
 from pwem.protocols import ProtProcessMovies
 from xmipp3.convert import getScipionObj
-import statistics as stat
-from pyworkflow import BETA, UPDATED, NEW, PROD
+from pyworkflow import NEW, PROD
 
 THRESHOLD = 2
 OUTPUT_ACCEPTED = 'outputMovies'
 OUTPUT_DISCARDED = 'outputMoviesDiscarded'
 
-class XmippProtMovieDoseAnalysis(ProtProcessMovies):
+class XmippProtMovieDoseAnalysis(ProtProcessMovies, Protocol):
     """ Protocol for the dose analysis """
     # FIXME: WITH .mrcs IT DOES NOT FILL THE LABELS
     _devStatus = PROD
-
     _label = 'movie dose analysis'
     _lastUpdateVersion = VERSION_3_0
-    _devStatus = NEW
     _possibleOutputs = {
         OUTPUT_ACCEPTED: SetOfMovies,
         OUTPUT_DISCARDED: SetOfMovies
@@ -311,24 +308,6 @@ def _checkNewOutput(self):
             if outputStep and outputStep.isWaiting():
                 outputStep.setStatus(cons.STATUS_NEW)
 
-    def _updateOutputSet(self, outputName, outputSet, state=Set.STREAM_OPEN):
-        outputSet.setStreamState(state)
-        if self.hasAttribute(outputName):
-            outputSet.write()  # Write to commit changes
-            outputAttr = getattr(self, outputName)
-            # Copy the properties to the object contained in the protocol
-            outputAttr.copy(outputSet, copyId=False)
-            # Persist changes
-            self._store(outputAttr)
-        else:
-            # Here the defineOutputs function will call the write() method
-            self._defineOutputs(**{outputName: outputSet})
-            self._store(outputSet)
-
-        # Close set databaset to avoid locking it
-        outputSet.close()
-
-
 # ------------------------- UTILS functions --------------------------------
     def getLimitIntervals(self):
         """ Funtion to obtain the acceptance interval limits."""
diff --git a/xmipp3/protocols/protocol_movie_opticalflow.py b/xmipp3/protocols/protocol_movie_opticalflow.py
index e6a07abbf..db34b242c 100644
--- a/xmipp3/protocols/protocol_movie_opticalflow.py
+++ b/xmipp3/protocols/protocol_movie_opticalflow.py
@@ -28,6 +28,7 @@
 # ******************************************************************************
 
 from os.path import join, exists
+import numpy as np
 
 import pyworkflow.object as pwobj
 import pyworkflow.utils as pwutils
@@ -245,7 +246,7 @@ def _validate(self):
             if self.doApplyDoseFilter:
                 doseFrame = inputSet.getAcquisition().getDosePerFrame()
 
-                if doseFrame == 0.0 or doseFrame is None:
+                if np.isclose(doseFrame, 0.0) or doseFrame is None:
                     errors.append('Dose per frame for input movies is 0 or not '
                                   'set. You cannot apply dose filter.')
 
diff --git a/xmipp3/protocols/protocol_multireference_alignability.py b/xmipp3/protocols/protocol_multireference_alignability.py
index ad6f615c7..adac5d97e 100644
--- a/xmipp3/protocols/protocol_multireference_alignability.py
+++ b/xmipp3/protocols/protocol_multireference_alignability.py
@@ -164,7 +164,7 @@ def _insertAllSteps(self):
             volDir = self._getVolDir(i + 1)
 
             pmStepId = self._insertFunctionStep('projectionLibraryStep',
-                                                volName, volDir,
+                                                volDir,
                                                 self.angularSampling.get(),
                                                 prerequisites=[convertId])
 
@@ -312,7 +312,7 @@ def phantomProject(self):
         self.runJob('xmipp_transform_mask', param, numberOfMpi=nproc,
                     numberOfThreads=nT)
 
-    def projectionLibraryStep(self, volName, volDir, angularSampling):
+    def projectionLibraryStep(self, volDir, angularSampling):
 
         # Generate projections from this reconstruction
         nproc = self.numberOfMpi.get()
@@ -574,10 +574,6 @@ def _defineMetadataRootName(self, mdrootname, volId):
             VolPrefix = 'vol%03d_' % (volId)
             return self._getExtraPath(VolPrefix + 'validation.xmd')
 
-    def _definePName(self):
-        fscFn = self._defineMetadataRootName('P')
-        return fscFn
-
     def _defineVolumeName(self, volId):
         fscFn = self._defineMetadataRootName('Volume', volId)
         return fscFn
diff --git a/xmipp3/protocols/protocol_particle_pick_automatic.py b/xmipp3/protocols/protocol_particle_pick_automatic.py
index 851b34a96..7f7e80aa4 100644
--- a/xmipp3/protocols/protocol_particle_pick_automatic.py
+++ b/xmipp3/protocols/protocol_particle_pick_automatic.py
@@ -47,7 +47,7 @@ class XmippParticlePickingAutomatic(ProtParticlePickingAuto, XmippProtocol):
     """Protocol to pick particles automatically in a set of
     micrographs using previous training """
     _label = 'auto-picking (step 2)'  
-    _devStatus = UPDATED
+    _devStatus = PROD
     filesToCopy = ['model_svm.txt', 'model_pca_model.stk', 'model_rotpca_model.stk',
                    'model_particle_avg.xmp', 'templates.stk']
     
diff --git a/xmipp3/protocols/protocol_particle_pick_consensus.py b/xmipp3/protocols/protocol_particle_pick_consensus.py
index f6e9df001..8992206b1 100644
--- a/xmipp3/protocols/protocol_particle_pick_consensus.py
+++ b/xmipp3/protocols/protocol_particle_pick_consensus.py
@@ -174,7 +174,8 @@ def _checkNewInput(self):
             # for non streaming do all and in the last iteration of streaming do the rest
             newMicIds = allMics.difference(self.checkedMics)
         else:  # for streaming processing, only go for the ready mics in all pickers
-            newMicIds = readyMics.difference(self.checkedMics)
+            if readyMics is not None:
+                newMicIds = readyMics.difference(self.checkedMics)
 
         if newMicIds:
             self.checkedMics.update(newMicIds)
diff --git a/xmipp3/protocols/protocol_particle_pick_remove_duplicates.py b/xmipp3/protocols/protocol_particle_pick_remove_duplicates.py
index ca0c9c045..165733808 100644
--- a/xmipp3/protocols/protocol_particle_pick_remove_duplicates.py
+++ b/xmipp3/protocols/protocol_particle_pick_remove_duplicates.py
@@ -36,6 +36,7 @@
 from .protocol_particle_pick_consensus import (XmippProtConsensusPicking,
                                                consensusWorker, getReadyMics)
 
+FACTOR_RADIUS = 0.6
 
 class XmippProtPickingRemoveDuplicates(XmippProtConsensusPicking):
     """
@@ -53,10 +54,12 @@ def _defineParams(self, form):
                       pointerClass='SetOfCoordinates',
                       label="Input coordinates", important=True,
                       help='Select the set of coordinates to compare')
-        form.addParam('consensusRadius', params.IntParam, default=10,
+        form.addParam('consensusRadius', params.IntParam, default=-1,
                       label="Radius",
+                      allowsPointers=True,
                       help="All coordinates within this radius (in pixels) "
-                           "are presumed to correspond to the same particle")
+                           "are presumed to correspond to the same particle.\n "
+                           "If -1 then 0.6 * box size is going to be assigned as radius.")
 
         # FIXME: It's not using more than one since
         #         self.stepsExecutionMode = STEPS_SERIAL
@@ -101,6 +104,13 @@ def _checkNewInput(self):
     def getMainInput(self):
         return self.inputCoordinates.get()
 
+    def getConsensusRadius(self):
+        if self.consensusRadius.get() == -1:
+            consensusRadius = int(self.inputCoordinates.get().getBoxSize()*FACTOR_RADIUS)
+        else:
+            consensusRadius = self.consensusRadius.get()
+        return consensusRadius
+
     def defineRelations(self, outputSet):
         self._defineTransformRelation(self.getMainInput(), outputSet)
 
@@ -112,7 +122,7 @@ def removeDuplicatesStep(self, micId, micName):
                                  self.getMainInput().iterCoordinates(micId)],
                                 dtype=int)
 
-        consensusWorker([coordArray], 1, self.consensusRadius.get(),
+        consensusWorker([coordArray], 1, self.getConsensusRadius(),
                         self._getTmpPath('%s%s.txt' % (self.FN_PREFIX, micId)))
 
         self.processedMics.update([micId])
@@ -122,4 +132,4 @@ def _validate(self):
         return errors
 
     def _summary(self):
-        return ["Radius = %d" % self.consensusRadius]
+        return ["Radius = %d" % self.getConsensusRadius()]
diff --git a/xmipp3/protocols/protocol_preprocess/protocol_process.py b/xmipp3/protocols/protocol_preprocess/protocol_process.py
index 4f83b1da6..07c949da0 100644
--- a/xmipp3/protocols/protocol_preprocess/protocol_process.py
+++ b/xmipp3/protocols/protocol_preprocess/protocol_process.py
@@ -157,7 +157,11 @@ def createOutputStep(self):
             vol.copyInfo(volInput)
             vol.setLocation(self.outputStk)
             if self.outputStk.endswith(".mrc"):
-                self.runJob("xmipp_image_header","-i %s --sampling_rate %f"%(self.outputStk,volInput.getSamplingRate()))
+                if hasattr(self, 'resizeOption') and self.doResize:
+                    samplingRate = self.resizeSamplingRate.get()
+                else:
+                    samplingRate = volInput.getSamplingRate()
+                self.runJob("xmipp_image_header","-i %s --sampling_rate %f"%(self.outputStk, samplingRate))
             if volInput.hasOrigin():
                 vol.setOrigin(volInput.getOrigin())
             self._postprocessOutput(vol)
diff --git a/xmipp3/protocols/protocol_projmatch/protocol_projmatch.py b/xmipp3/protocols/protocol_projmatch/protocol_projmatch.py
index d702a2919..ec8a50e18 100644
--- a/xmipp3/protocols/protocol_projmatch/protocol_projmatch.py
+++ b/xmipp3/protocols/protocol_projmatch/protocol_projmatch.py
@@ -27,8 +27,8 @@
 # **************************************************************************
 
 from pyworkflow.object import Integer
-from pyworkflow.utils import (getFloatListFromValues, getBoolListFromValues,
-                              getStringListFromValues)
+from pyworkflow.utils import (getBoolListFromValues,
+                              getStringListFromValues, getListFromValues)
 from pwem.protocols import ProtRefine3D, ProtClassify3D
 from pwem import emlib
 from .projmatch_initialize import *
@@ -267,7 +267,7 @@ def itersFloatValues(self, attributeName, firstValue=-1):
         valuesStr = self.getAttributeValue(attributeName)
         if valuesStr is None:
             raise Exception('None value for attribute: %s' % attributeName)
-        return [firstValue] + getFloatListFromValues(valuesStr, length=self.numberOfIterations.get())
+        return [firstValue] + getListFromValues(valuesStr, length=self.numberOfIterations.get(), caster=float)
     
     def itersBoolValues(self, attributeName, firstValue=False):
         """ Take the string of a given attribute and
diff --git a/xmipp3/protocols/protocol_reconstruct_fourier.py b/xmipp3/protocols/protocol_reconstruct_fourier.py
index a9125b193..e80641f86 100644
--- a/xmipp3/protocols/protocol_reconstruct_fourier.py
+++ b/xmipp3/protocols/protocol_reconstruct_fourier.py
@@ -26,12 +26,14 @@
 
 from pwem.objects import Volume
 from pwem.protocols import ProtReconstruct3D
+from pwem import emlib
 import pyworkflow.protocol.params as params
 import pyworkflow.protocol.constants as cons
 from xmipp3.convert import writeSetOfParticles
 from xmipp3.base import isXmippCudaPresent
+from pyworkflow.utils import moveFile
 import os
-
+from pyworkflow import BETA, UPDATED, NEW, PROD
 
 class XmippProtReconstructFourier(ProtReconstruct3D):
     """    
@@ -40,7 +42,8 @@ class XmippProtReconstructFourier(ProtReconstruct3D):
     and used as direction projections to reconstruct.
     """
     _label = 'reconstruct fourier'
-    
+    _devStatus = UPDATED
+
     #--------------------------- DEFINE param functions --------------------------------------------   
     def _defineParams(self, form):
 
@@ -68,6 +71,8 @@ def _defineParams(self, form):
                       help='Maximum resolution (in Angstrom) to consider \n'
                            'in Fourier space (default Nyquist).\n'
                            'Param *--maxres* in Xmipp.')
+        form.addParam('useHalves', params.BooleanParam, label='Use halves', default=False,
+                      help='Create separate reconstructions from two random subsets. Useful for resolution measurements')
         line = form.addLine('Padding factor',
                              expertLevel=cons.LEVEL_ADVANCED,
                              help='Padding of the input images. Higher number will result in more precise interpolation in Fourier '
@@ -102,21 +107,36 @@ def _createFilenameTemplates(self):
         """ Centralize how files are called for iterations and references. """
         myDict = {
             'input_xmd': self._getExtraPath('input_particles.xmd'),
-            'output_volume': self._getPath('output_volume.mrc')
+            'half1_xmd': self._getExtraPath('input_particles000001.xmd'),
+            'half2_xmd': self._getExtraPath('input_particles000002.xmd'),
+            'output_volume': self._getPath('output_volume.mrc'),
+            'half1_volume': self._getPath('half1.mrc'),
+            'half2_volume': self._getPath('half2.mrc')
             }
         self._updateFilenamesDict(myDict)
 
     def _insertAllSteps(self):
         self._createFilenameTemplates()
         self._insertFunctionStep('convertInputStep')
-        self._insertReconstructStep()
+        if self.useHalves.get():
+            self._insertFunctionStep('splitInputStep')
+            self._insertReconstructStep('half1')
+            self._insertReconstructStep('half2')
+            self._insertFunctionStep('averageStep')
+        else:
+            self._insertReconstructStep()
         self._insertFunctionStep('createOutputStep')
         
-    def _insertReconstructStep(self):
+    def _insertReconstructStep(self, half=None):
         #imgSet = self.inputParticles.get()
 
-        params =  '  -i %s' % self._getFileName('input_xmd')
-        params += '  -o %s' % self._getFileName('output_volume')
+        if half is None:
+            params =  '  -i %s' % self._getFileName('input_xmd')
+            params += '  -o %s' % self._getFileName('output_volume')
+        else:
+            params =  '  -i %s' % self._getFileName(half + '_xmd')
+            params += '  -o %s' % self._getFileName(half + '_volume')
+            
         params += ' --sym %s' % self.symmetryGroup.get()
         maxRes = self.maxRes.get()
         if maxRes == -1:
@@ -163,7 +183,14 @@ def convertInputStep(self):
         #TODO: This only writes metadata what about binary file
         #it should
         writeSetOfParticles(imgSet, particlesMd)
+        
+    def splitInputStep(self):
+        args = []
+        args += ['-i', self._getFileName('input_xmd')]
+        args += ['-n', 2]
 
+        self.runJob('xmipp_metadata_split', args, numberOfMpi=1)
+        
     def reconstructStep(self, params):
         """ Create the input file in STAR format as expected by Xmipp.
         If the input particles comes from Xmipp, just link the file. 
@@ -178,16 +205,34 @@ def reconstructStep(self, params):
                 self.runJob('xmipp_reconstruct_fourier', params)
             else:
                 self.runJob('xmipp_reconstruct_fourier_accel', params)
-
-        self.runJob("xmipp_image_header", "-i %s --sampling_rate %f"%\
-                    (self._getFileName('output_volume'), self.inputParticles.get().getSamplingRate()),
-                    numberOfMpi=1)
             
+    def averageStep(self):
+        # Read
+        half1 = emlib.Image(self._getFileName('half1_volume'))
+        half2 = emlib.Image(self._getFileName('half2_volume'))
+        
+        # Average
+        half1.inplaceAdd(half2)
+        half1.inplaceMultiply(0.5)
+        
+        # Write
+        half1.write(self._getFileName('output_volume'))
+        
     def createOutputStep(self):
         imgSet = self.inputParticles.get()
+
+        self.runJob("xmipp_image_header", "-i %s --sampling_rate %f"%\
+                    (self._getFileName('output_volume'), imgSet.getSamplingRate()),
+                    numberOfMpi=1)
+        
         volume = Volume()
         volume.setFileName(self._getFileName('output_volume'))
         volume.setSamplingRate(imgSet.getSamplingRate())
+        if self.useHalves.get():
+            volume.setHalfMaps([
+                self._getFileName('half1_volume'),
+                self._getFileName('half2_volume')
+            ])
         
         self._defineOutputs(outputVolume=volume)
         self._defineSourceRelation(self.inputParticles, volume)
diff --git a/xmipp3/protocols/protocol_resolution_bfactor.py b/xmipp3/protocols/protocol_resolution_bfactor.py
index 786bea3eb..cd8bcd0ca 100644
--- a/xmipp3/protocols/protocol_resolution_bfactor.py
+++ b/xmipp3/protocols/protocol_resolution_bfactor.py
@@ -41,7 +41,7 @@
                                         LEVEL_ADVANCED)
 from pyworkflow.object import Float
 from pyworkflow.utils import getExt
-from pwem.objects import Volume, SetOfParticles
+from pwem.objects import Volume, SetOfParticles, AtomStruct
 
 
 FN_METADATA_BFACTOR_RESOLUTION = 'bfactor_resolution.xmd'
@@ -70,13 +70,13 @@ def _defineParams(self, form):
 
         form.addParam('normalizeResolution', BooleanParam, default=True,
                       label="Normalize Resolution",
-                      help='The normalizedlocal resolution map is defined as'
-                           '(LR - FSC)/FSC, where LR is the local resolution of a'
+                      help='The normalizedlocal resolution map is defined as '
+                           '(LR - FSC)/FSC, where LR is the local resolution of a '
                            'given voxel, and FSC is the FSC resolution in A. This '
-                           'map provides information about whether the local resolution is'
+                           'map provides information about whether the local resolution is '
                            'greater or lesser than the FSC. The local resolution '
-                           'normalized map is used to carry out the matching with the local'
-                           'bfactor per residue. Yes means that the local resolution will be'
+                           'normalized map is used to carry out the matching with the local '
+                           'bfactor per residue. Yes means that the local resolution will be '
                            'normalized by the algorithm. No means that the input local '
                            'resolution map is already a normalized local resolution map.')
 
@@ -85,7 +85,7 @@ def _defineParams(self, form):
                       condition='normalizeResolution',
                       help='Select a local resolution map. Alternatively, the input.'
                            ' can be a normalized local resolution map, in this case'
-                           'set the Normalize resolution to No')
+                           ' set the Normalize resolution to No')
 
         form.addParam('normalizedMap', PointerParam, pointerClass='Volume',
                       label="Normalized Local Resolution Map", important=True,
@@ -93,7 +93,7 @@ def _defineParams(self, form):
                       help='Select a normalized local resolution map. The local'
                            ' resolution normalized map is defined as '
                            '(LR - FSC)/FSC, where LR is the local resolution of a'
-                           'given voxel, and FSC is the FSC resolution in A')
+                           ' given voxel, and FSC is the FSC resolution in A')
 
         form.addParam('fscResolution', FloatParam,
                       condition = 'normalizeResolution',
@@ -102,7 +102,7 @@ def _defineParams(self, form):
 
         form.addParam('medianEstimation', BooleanParam, default=True,
                       label="Use median",
-                      help='The local resolution per residue can be estimated using'
+                      help='The local resolution per residue can be estimated using '
                            'the mean (by default - No) or the median (yes)')
 
         form.addParam('centered', BooleanParam, default=True,
@@ -122,6 +122,8 @@ def _insertAllSteps(self):
         # 2 Carry out the matching betweent the local resolution per residue and bfactor
         self._insertFunctionStep('matchingBfactorLocalResolution')
 
+        self._insertFunctionStep('createOutputStep')
+
     def mrc_convert(self, fileName, outputFileName):
         """Check if the extension is .mrc, if not then uses xmipp to convert it
         """
@@ -137,8 +139,13 @@ def mrc_convert(self, fileName, outputFileName):
     def convertInputStep(self):
         """ Read the input volume and check the file extension to convert to mrc is it is the case.
         """
-        self.vol = self.mrc_convert(self.localResolutionMap.get().getFileName(),
-                                    self._getTmpPath('localResolutionMap.mrc'))
+        if self.normalizeResolution.get():
+            self.inputMap = self.localResolutionMap
+        else:
+            self.inputMap = self.normalizedMap
+
+        self.vol = self.mrc_convert(self.inputMap.get().getFileName(),
+                                    self._getTmpPath('localResolutionMap.mrc'))        
 
     def matchingBfactorLocalResolution(self):
         """ The local resolution map and the pdb are taken and analyzed to match the
@@ -151,7 +158,7 @@ def matchingBfactorLocalResolution(self):
         params += ' --vol %s' % self.vol
         if self.normalizeResolution.get():
             params += ' --fscResolution %s' % self.fscResolution.get()
-        params += ' --sampling %f' % self.localResolutionMap.get().getSamplingRate()
+        params += ' --sampling %f' % self.inputMap.get().getSamplingRate()
         if self.medianEstimation.get():
             params += ' --useMedian '
         if self.centered.get():
@@ -159,3 +166,8 @@ def matchingBfactorLocalResolution(self):
         params += ' -o %s' % self._getExtraPath()
 
         self.runJob('xmipp_resolution_pdb_bfactor', params)
+
+    def createOutputStep(self):
+        outputPdb = AtomStruct()
+        outputPdb.setFileName(self._getExtraPath("chimeraPDB.pdb"))
+        self._defineOutputs(outputStructure=outputPdb)
diff --git a/xmipp3/protocols/protocol_resolution_deepres.py b/xmipp3/protocols/protocol_resolution_deepres.py
index b2064434e..9ec2f7268 100644
--- a/xmipp3/protocols/protocol_resolution_deepres.py
+++ b/xmipp3/protocols/protocol_resolution_deepres.py
@@ -124,8 +124,9 @@ def _createFilenameTemplates(self):
 
     def _insertAllSteps(self):
             # Convert input into xmipp Metadata format
-        
-        updateEnviron(self.gpuList.get())    
+
+        if not self.useQueueForSteps() and not self.useQueue():
+            updateEnviron(self.gpuList.get())
             
         self._createFilenameTemplates() 
         self._insertFunctionStep('convertInputStep')
diff --git a/xmipp3/protocols/protocol_screen_deepConsensus.py b/xmipp3/protocols/protocol_screen_deepConsensus.py
index 99e74bbdf..7bac83a3a 100644
--- a/xmipp3/protocols/protocol_screen_deepConsensus.py
+++ b/xmipp3/protocols/protocol_screen_deepConsensus.py
@@ -1106,9 +1106,9 @@ def predictCNN(self):
 
         args= " -n %s --mode score -i %s -o %s "%(netDataPath, fnamesPred, outParticlesPath)
 
-        if posTestDict and posTestDict:
+        if posTestDict and negTestDict:
           fnamesPosTest, weightsPosTest= self.__dataDict_toStrs(posTestDict)
-          fnamesNegTest, weightsNegTest= self.__dataDict_toStrs(posTestDict)
+          fnamesNegTest, weightsNegTest= self.__dataDict_toStrs(negTestDict)
           args+= " --testingTrue %s --testingFalse %s "%(fnamesPosTest, fnamesNegTest)
 
         if not gpuToUse is None:
diff --git a/xmipp3/protocols/protocol_screen_particles.py b/xmipp3/protocols/protocol_screen_particles.py
index 9a86ca5a5..585605b96 100644
--- a/xmipp3/protocols/protocol_screen_particles.py
+++ b/xmipp3/protocols/protocol_screen_particles.py
@@ -364,7 +364,8 @@ def _summary(self):
             summary += sumRejMet.values()
             summary.append("Output particles not ready yet.")
         else:
-            summary.append(sumRejMet['Zscore'])
+            if 'Zscore' in sumRejMet:
+                summary.append(sumRejMet['Zscore'])
             if hasattr(self, 'sumZScore'):
                 summary.append(" - The minimum ZScore is %.2f" % self.minZScore)
                 summary.append(" - The maximum ZScore is %.2f" % self.maxZScore)
@@ -373,8 +374,10 @@ def _summary(self):
             else:
                 summary.append(
                     "Summary values not calculated during processing.")
-            summary.append(sumRejMet['SSNR'])
-            summary.append(sumRejMet['Var'])
+            if 'SSNR' in sumRejMet:
+                summary.append(sumRejMet['SSNR'])
+            if 'Var' in sumRejMet:
+                summary.append(sumRejMet['Var'])
             if self.autoParRejectionVar != self.REJ_NONE:
                 if hasattr(self, 'varThreshold'):
                     summary.append(" - Variance threshold: %.2f" % self.varThreshold)
diff --git a/xmipp3/protocols/protocol_simulate_ctf.py b/xmipp3/protocols/protocol_simulate_ctf.py
index 7ea1fb694..1651d9dc6 100644
--- a/xmipp3/protocols/protocol_simulate_ctf.py
+++ b/xmipp3/protocols/protocol_simulate_ctf.py
@@ -25,14 +25,16 @@
 # *
 # **************************************************************************
 
+import numpy as np
 import random
 
 import pyworkflow.protocol.params as params
 from pyworkflow import VERSION_3_0
 from pwem.protocols import Prot2D
 from pwem.objects import CTFModel
+from pyworkflow import BETA, UPDATED, NEW, PROD
 
-from xmippLib import createEmptyFile
+import xmippLib
 
 
 class XmippProtSimulateCTF(Prot2D):
@@ -42,6 +44,8 @@ class XmippProtSimulateCTF(Prot2D):
     """
 
     _label = 'simulate ctf'
+    _devStatus = UPDATED
+
     _lastUpdateVersion = VERSION_3_0
 
     def __init__(self, *args, **kwargs):
@@ -81,6 +85,8 @@ def _defineParams(self, form):
                       label="Lower defocus difference between defocusU and defocusV (A)")
         form.addParam('DefocusFdiff', params.FloatParam, default=500, condition='astig',
                       label="Upper defocus difference between defocusU and defocusV(A)")
+        form.addParam('noiseBefore', params.FloatParam, default=0, label='Noise before CTF', help='Sigma')
+        form.addParam('noiseAfter', params.FloatParam, default=0, label='Noise after CTF', help='Sigma')
 
     # --------------------------- INSERT steps functions ------------------------
     def _insertAllSteps(self):
@@ -91,7 +97,7 @@ def _insertAllSteps(self):
     def convertInputStep(self):
         x, y, _ = self.inputParticles.get().getDimensions()
         n = self.inputParticles.get().getSize()
-        createEmptyFile(self._getPath("images.mrc"), x, y, 1, n)
+        xmippLib.createEmptyFile(self._getPath("images.mrc"), x, y, 1, n)
 
     def simulateStep(self):
         n = 1
@@ -111,6 +117,14 @@ def simulateStep(self):
             location = particle.getLocation()
             fnIn = str(location[0]) + "@" + location[1]
             fnOut = str(n) + "@" + fnStk
+
+            if self.noiseBefore>0:
+                I=xmippLib.Image(fnIn)
+                Idata = I.getData()
+                I.setData(Idata+self.noiseBefore.get()*np.numpy.random.Generator.normal(size=Idata.shape))
+                I.write(fnOut)
+                fnIn=fnOut
+
             defocusU = random.uniform(self.Defocus0.get(), self.DefocusF.get())
             args = "-i %s -o %s" % (fnIn, fnOut)
             if self.astig:
@@ -125,6 +139,12 @@ def simulateStep(self):
                         (self.voltage, self.cs, self.Q0, defocusU, Ts)
             self.runJob("xmipp_transform_filter", args)
 
+            if self.noiseAfter>0:
+                I=xmippLib.Image(fnOut)
+                Idata = I.getData()
+                I.setData(Idata+self.noiseAfter.get()*np.random.Generator.normal(size=Idata.shape))
+                I.write(fnOut)
+
             newCTF = CTFModel()
             newCTF.setDefocusU(defocusU)
             newCTF.setDefocusV(defocusV)
@@ -150,5 +170,7 @@ def _summary(self):
         summary.append("Cs=%f mm" % self.cs)
         summary.append("Q0=%f" % self.Q0)
         summary.append("Defocus range=[%f,%f] A" % (self.Defocus0, self.DefocusF))
+        summary.append('Noise before=%f' %self.noiseBefore)
+        summary.append('Noise after=%f' %self.noiseAfter)
         return summary
 
diff --git a/xmipp3/protocols/protocol_tilt_analysis.py b/xmipp3/protocols/protocol_tilt_analysis.py
index debb95d27..a4b13504e 100644
--- a/xmipp3/protocols/protocol_tilt_analysis.py
+++ b/xmipp3/protocols/protocol_tilt_analysis.py
@@ -44,8 +44,8 @@
 from pwem.protocols import ProtMicrographs
 from xmipp3 import emlib
 from xmipp3.convert import getScipionObj
-from matplotlib import pyplot as plt
 
+AUTOMATIC_WINDOW_SIZES = [4096, 2048, 1024, 512, 256]
 
 class XmippProtTiltAnalysis(ProtMicrographs, Protocol):
     """ Estimate the tilt of a micrograph, by analyzing the PSD correlations of different segments of the image.
@@ -67,10 +67,17 @@ def _defineParams(self, form):
                       label="Input micrographs", important=True,
                       help='Select the SetOfMicrograph to be preprocessed.')
 
+        form.addParam('autoWindow', BooleanParam, default=True, expertLevel=LEVEL_ADVANCED,
+                      label='Estimate automatically the window size?',
+                      help='Use this button to decide if you want to estimate the window size'
+                           'based on the movie size. It will select between 512, 1024, 2048 or 4096. '
+                           'Condition:  window_size <= Size/2.5 ')
+
         form.addParam('window_size', IntParam, label='Window size',
+                      condition="not autoWindow",
                       default=1024, expertLevel=LEVEL_ADVANCED, validators=[GE(256, 'Error must be greater than 256')],
-                      help='''By default, the micrograph will be divided into windows of dimensions 512x512, the PSD '''
-                           '''its correlations will be computed in every segment.''')
+                      help=''' By default, the micrograph will be divided into windows of dimensions 1024x1024, '''
+                           ''' the PSD its correlations will be computed in every segment.''')
 
         form.addParam('objective_resolution', FloatParam, label='Objective resolution',
                       default=3, expertLevel=LEVEL_ADVANCED, validators=[GT(0, 'Error must be Positive')],
@@ -118,6 +125,7 @@ def initializeStep(self):
         self.streamClosed = self.inputMicrographs.get().isStreamClosed()
         self.micsDict = {mic.getObjId(): mic.clone() for mic in self._loadInputList(self.micsFn).iterItems()}
         pwutils.makePath(self._getExtraPath('DONE'))
+        self.windowSize = self.getWindowSize()
 
     def createOutputStep(self):
         self._closeOutputSet()
@@ -337,9 +345,9 @@ def calculateTiltCorrelationStep(self, mic):
         micFolder = self._getOutputMicFolder(mic)
         micImage = ImageHandler().read(mic.getLocation())
         dimx, dimy, z, n = micImage.getDimensions()
-        wind_step = self.window_size.get()
-        x_steps, y_steps = window_coordinates2D(dimx, dimy, wind_step)
-        subWindStep = int(wind_step * (self.samplingRate / self.objective_resolution.get()))
+        windStep = self.windowSize
+        x_steps, y_steps = window_coordinates2D(dimx, dimy, windStep)
+        subWindStep = int(windStep * (self.samplingRate / self.objective_resolution.get()))
         x_steps_psd, y_steps_psd = window_coordinates2D(subWindStep * 2, subWindStep * 2, subWindStep)
         # Extract windows
         window_image = ImageHandler().createImage()
@@ -472,6 +480,23 @@ def _correctFormat(self, micName, micFn, micFolderTmp):
 
         return newMicName
 
+    def getWindowSize(self):
+        """ Function to get the window size, automatically or the one set by the user. """
+        if self.autoWindow:
+            dimX, dimY, _ = self.getInputMicrographs().getFirstItem().getDimensions()
+            halfMic = int(min(dimX, dimY)/2.5)
+            windowSize = halfMic  # In case there is not a suitable option, very unlikely
+            for sizeAuto in AUTOMATIC_WINDOW_SIZES:
+                if sizeAuto < halfMic:  # Exit in case you found a suitable option
+                    windowSize = sizeAuto
+                    break
+        else:
+            windowSize = self.window_size.get()
+
+        self.info('The window size used is %d' % windowSize)
+
+        return windowSize
+
     def _insertFinalSteps(self, deps):
         """ This should be implemented in subclasses"""
         return deps
@@ -568,7 +593,7 @@ def rotation(imag, angle, shape, P):
     return transformed, M
 
 
-def window_coordinates2D(x, y, wind_step):
+def window_coordinates2D(x, y, windStep):
     x0 = 0
     xF = x - 1
     y0 = 0
@@ -576,15 +601,15 @@ def window_coordinates2D(x, y, wind_step):
     x_coor = []
     y_coor = []
 
-    if wind_step < x and wind_step < y:
+    if windStep < x and windStep < y:
         x_coor.append(x0)
-        x_coor.append(x0 + wind_step - 1)
-        x_coor.append(xF - wind_step)
+        x_coor.append(x0 + windStep - 1)
+        x_coor.append(xF - windStep)
         x_coor.append(xF)
 
         y_coor.append(y0)
-        y_coor.append(y0 + wind_step - 1)
-        y_coor.append(yF - wind_step)
+        y_coor.append(y0 + windStep - 1)
+        y_coor.append(yF - windStep)
         y_coor.append(yF)
 
         return x_coor, y_coor
diff --git a/xmipp3/protocols/protocol_validate_fscq.py b/xmipp3/protocols/protocol_validate_fscq.py
index 9a6e9bb20..af058885f 100644
--- a/xmipp3/protocols/protocol_validate_fscq.py
+++ b/xmipp3/protocols/protocol_validate_fscq.py
@@ -25,7 +25,7 @@
 # *
 # **************************************************************************
 
-import os
+import os, shutil, gzip
 from pyworkflow import VERSION_3_0
 from pyworkflow.protocol.params import (PointerParam, BooleanParam,
                                         IntParam, FileParam, FloatParam)
@@ -38,6 +38,7 @@
 from pwem.convert import Ccp4Header
 from pwem.convert.atom_struct import toPdb, toCIF, AtomicStructHandler, addScipionAttribute
 import xmipp3
+from pyworkflow import BETA, UPDATED, NEW, PROD
 
 
 VALIDATE_METHOD_URL = 'https://github.com/I2PC/scipion-em-xmipp/wiki/XmippProtValFit'
@@ -47,7 +48,6 @@
 BLOCRES_HALF_FILE = 'blocresHalf'
 RESTA_FILE = 'diferencia.vol'
 RESTA_FILE_MRC = 'diferencia.map'
-PDB_VALUE_FILE = 'pdb_fsc-q.pdb'
 MASK_FILE_MRC = 'mask.map'
 MASK_FILE = 'mask.vol'
 FN_VOL = 'vol.map'
@@ -56,7 +56,6 @@
 MD_MEANS = 'params.xmd'
 MD2_MEANS = 'params2.xmd'
 RESTA_FILE_NORM = 'diferencia_norm.map'
-PDB_NORM_FILE = 'pdb_fsc-q_norm.pdb'
 OUTPUT_CIF = 'fscq_struct.cif'
 
 
@@ -66,6 +65,7 @@ class XmippProtValFit(ProtAnalysis3D):
     """
     _label = 'validate fsc-q'
     _lastUpdateVersion = VERSION_3_0
+    _devStatus = UPDATED
     _ATTRNAME = 'fscq_score'
     _OUTNAME = 'outputAtomStruct'
     _possibleOutputs = {_OUTNAME: AtomStruct}
@@ -127,7 +127,6 @@ def _defineParams(self, form):
         form.addParallelSection(threads=8, mpi=1)
 
     # --------------------------- INSERT steps functions --------------------------------------------
-
     def _createFilenameTemplates(self):
         """ Centralize how files are called """
         myDict = {
@@ -138,8 +137,6 @@ def _createFilenameTemplates(self):
                  RESTA_FILE: self._getTmpPath('diferencia.vol'),
                  RESTA_FILE_MRC: self._getExtraPath('diferencia.map'),
                  RESTA_FILE_NORM: self._getExtraPath('diferencia_norm.map'),
-                 PDB_VALUE_FILE:  self._getExtraPath('pdb_fsc-q.pdb'),
-                 PDB_NORM_FILE: self._getExtraPath('pdb_fsc-q_norm.pdb'),
                  MASK_FILE_MRC : self._getExtraPath('mask.map'),
                  MASK_FILE: self._getTmpPath('mask.vol'),
                  FN_VOL: self._getTmpPath("vol.map"),
@@ -165,10 +162,10 @@ def _insertAllSteps(self):
     def convertInputStep(self):
         """ Convert inputs to desired format."""
         #Convert Input to pdb
-        if self.getInputStructFile().endswith('.pdb') or self.getInputStructFile().endswith('.ent'):
-            os.symlink(os.path.abspath(self.getInputStructFile()), self.getPDBFile())
+        if self.isStructExtensionValid():
+            os.symlink(os.path.abspath(self.getInputStructFile()), self.getStructFile())
         else:
-            toPdb(self.getInputStructFile(), self.getPDBFile())
+            toPdb(self.getInputStructFile(), self.getStructFile())
 
         """ Read the input volume."""
         self.volume = self.inputVolume.get()
@@ -221,7 +218,7 @@ def convertInputStep(self):
             params += ' -v 0 '
             params += ' --sampling %f' % self.inputVolume.get().getSamplingRate()
             params += ' --size %d' % self.inputVolume.get().getXDim()
-            params += ' -i %s' % self.getPDBFile()
+            params += ' -i %s' % self.getStructFile()
             params += ' -o %s' % self._getFileName(OUTPUT_PDBVOL_FILE)
             self.runJob('xmipp_volume_from_pdb', params)
 
@@ -363,10 +360,10 @@ def substractBlocresStep(self):
 
     def assignPdbStep(self):
 
-        params = ' --pdb %s ' % self.getPDBFile()
+        params = ' --pdb %s ' % self.getStructFile()
         params += ' --vol %s ' % self._getFileName(RESTA_FILE_MRC)
         params += ' --mask %s ' % self.mask_xmipp
-        params += ' -o %s ' % self._getFileName(PDB_VALUE_FILE)
+        params += ' -o %s ' % self.getFSCQFile()
         params += ' --sampling %f' % self.inputVolume.get().getSamplingRate()
         params += ' --origin %f %f %f' %(self.x, self.y, self.z)
         params += ' --radius 0.8'
@@ -374,10 +371,10 @@ def assignPdbStep(self):
         self.runJob('xmipp_pdb_label_from_volume', params)
 
         """Diveded by resolution"""
-        params = ' --pdb %s ' % self.getPDBFile()
+        params = ' --pdb %s ' % self.getStructFile()
         params += ' --vol %s ' % self._getFileName(RESTA_FILE_NORM)
         params += ' --mask %s ' % self.mask_xmipp
-        params += ' -o %s ' % self._getFileName(PDB_NORM_FILE)
+        params += ' -o %s ' % self.getNormFSCQFile()
         params += ' --sampling %f' % self.inputVolume.get().getSamplingRate()
         params += ' --origin %f %f %f' %(self.x, self.y, self.z)
         params += ' --radius 0.8'
@@ -403,48 +400,55 @@ def _setMetrics(self):
 
         # statistic from fnal pdb with fsc-q
         # Number of atoms greater or less than 0.5
-        total_atom = 0
-        fscq_greater = 0
-        fscq_less = 0
-        with open(self._getFileName(PDB_VALUE_FILE)) as f:
-            lines_data = f.readlines()
-            for j, lin in enumerate(lines_data):
+        totalAtom = 0
+        fscqGreater = 0
+        fscqLess = 0
 
-                if (lin.startswith('ATOM') or lin.startswith('HETATM')):
+        # Reading FCSQ's value file
+        atStHandler = AtomicStructHandler()
+        atStHandler.read(self.getFSCQFile())
 
-                    total_atom = total_atom + 1
-                    fscq_atom = float(lin[54:60])
+        # Reading value stored in occupancy field for each atom
+        for model in atStHandler.structure:
+            for atom in model.get_atoms():
+                totalAtom += 1
+                fscqAtom = atom.get_occupancy()
 
-                    if (fscq_atom > 0.5):
-                        fscq_greater = fscq_greater + 1
+                if (fscqAtom > 0.5):
+                    fscqGreater += 1
 
-                    if (fscq_atom < -0.5):
-                        fscq_less = fscq_less + 1
+                if (fscqAtom < -0.5):
+                    fscqLess += 1
 
-        porc_greater = (fscq_greater * 100) / total_atom
-        porc_less = (fscq_less * 100) / total_atom
+        porcGreater = (fscqGreater * 100) / totalAtom
+        porcLess = (fscqLess * 100) / totalAtom
 
         self.mean = Float(mean)
         self.meanA = Float(meanA)
         self.mean2 = Float(mean2)
         self.meanA2 = Float(meanA2)
-        self.total_atom = Integer(total_atom)
-        self.fscq_greater = Integer(fscq_greater)
-        self.fscq_less = Integer(fscq_less)
-        self.porc_greater = Float(porc_greater)
-        self.porc_less = Float(porc_less)
+        self.total_atom = Integer(totalAtom)
+        self.fscq_greater = Integer(fscqGreater)
+        self.fscq_less = Integer(fscqLess)
+        self.porc_greater = Float(porcGreater)
+        self.porc_less = Float(porcLess)
         self._store()
 
     def getFscqAttrDic(self):
         fscqDic = {}
-        with open(self._getFileName(PDB_VALUE_FILE)) as f:
-            lines_data = f.readlines()
-            for j, lin in enumerate(lines_data):
-                if (lin.startswith('ATOM') or lin.startswith('HETATM')):
-                    resNumber, resChain, atomName = lin[22:26].strip(), lin[21].strip(), lin[12:16].strip()
-                    resId = '{}:{}@{}'.format(resChain, resNumber, atomName)
-                    fscq_atom = lin[54:60].strip()
-                    fscqDic[resId] = fscq_atom
+
+        # Reading FCSQ's value file
+        atStHandler = AtomicStructHandler()
+        atStHandler.read(self.getFSCQFile())
+
+        # Reading value stored in occupancy field for each atom
+        for model in atStHandler.structure:
+            for atom in model.get_atoms():
+                fId = atom.get_full_id()
+                chainName, resNumber, atomName = fId[2], fId[3][1], fId[4][0]
+                resId = '{}:{}@{}'.format(chainName, resNumber, atomName)
+                fscqDic[resId] = atom.get_occupancy() # FCSQ value (stored in occupancy variable)
+
         return fscqDic
 
     def createOutputStep(self):
@@ -452,9 +456,9 @@ def createOutputStep(self):
         fscqDic = self.getFscqAttrDic()
 
         AS = self.getInputStruct()
+        ASCIF = self.getInputStructCIF()
         ASH = AtomicStructHandler()
-        inpAS = toCIF(AS.getFileName(), self._getTmpPath('inputStruct.cif'))
-        cifDic = ASH.readLowLevel(inpAS)
+        cifDic = ASH.readLowLevel(ASCIF)
         cifDic = addScipionAttribute(cifDic, fscqDic, self._ATTRNAME, recipient='atoms')
         ASH._writeLowLevel(self._getPath(OUTPUT_CIF), cifDic)
 
@@ -538,9 +542,40 @@ def getInputStructFile(self):
             return self.inputPDB.get()
         else:
             return self.inputPDBObj.get().getFileName()
-
-    def getPDBFile(self):
-        return self._getExtraPath('inputPDB.pdb')
+    
+    def getInputStructCIF(self) -> str:
+        """
+        ### This function returns the full path for the input Atom Struct file in cif format.
+        #### If file is not in cif format, it is converted to it. 
+
+        #### Returns:
+        - (str): Input Atom Struct's file's path.
+        """
+        # Get raw input filename
+        structFile = self.getInputStructFile()
+
+        if '.cif' not in structFile:
+            # If file is not in cif format, convert it
+            structFile = toCIF(structFile, self._getTmpPath('inputStruct.cif'))
+        elif structFile.endswith('.cif.gz'):
+            # If file is cif but compressed as a gz, extract it
+            oFile = self._getTmpPath('inputStruct.cif')
+            with gzip.open(structFile, 'rb') as fIn:
+                with open(oFile, 'wb') as fOut:
+                    shutil.copyfileobj(fIn, fOut)
+            structFile = oFile
+        
+        # Return input cif file
+        return structFile
+
+    def getStructFile(self) -> str:
+        """
+        ### This function returns the full path for the Atom Struct file.
+
+        #### Returns:
+        - (str): Atom Struct's file's path.
+        """
+        return self._getExtraPath('inputStruct' + self.getStructExtension())
 
     def getInputStruct(self):
         if self.fromFile:
@@ -548,3 +583,51 @@ def getInputStruct(self):
         else:
             return self.inputPDBObj.get()
 
+    def isStructExtensionValid(self) -> bool:
+        """
+        ### This function returns True if the input struct extension is one of the accepted types.
+
+        #### Returns:
+        - (bool): True if the input struct has an accepted extension. False otherwise.
+        """
+        # Getting file's natural extension
+        extension = os.path.splitext(self.getInputStructFile())[-1]
+
+        # Return extension validity
+        return extension == '.pdb' or extension == '.ent' or extension == '.cif' or extension == ".mmcif" or extension == ".cif.gz"
+    
+    def getStructExtension(self) -> str:
+        """
+        ### This function returns the expected extension for the Atom Struct file.
+
+        #### Returns:
+        - (str): Atom Struct's file's extension.
+        """
+        # Getting file's natural extension
+        extension = os.path.splitext(self.getInputStructFile())[-1]
+
+        # If extension is a maintainable type, return as is (not all valid types are maintainable)
+        # .pdb is not considered "maintainable" type as it is the default
+        if extension == '.cif' or extension == ".mmcif" or extension == ".cif.gz":
+            return extension
+        
+        # Return .pdb by default
+        return '.pdb'
+    
+    def getFSCQFile(self) -> str:
+        """
+        ### This function returns the filename for FSCQ's value file.
+
+        #### Returns:
+        - (str): FSCQ's value file's full name with extension.
+        """
+        return self._getExtraPath(f'pdb_fsc-q{self.getStructExtension()}')
+
+    def getNormFSCQFile(self) -> str:
+        """
+        ### This function returns the filename for FSCQ's norm file.
+
+        #### Returns:
+        - (str): FSCQ's norm file's full name with extension.
+        """
+        return self._getExtraPath(f'pdb_fsc-q_norm{self.getStructExtension()}')
diff --git a/xmipp3/protocols/protocol_validate_nontilt.py b/xmipp3/protocols/protocol_validate_nontilt.py
index 31de90eb8..0365df7c8 100644
--- a/xmipp3/protocols/protocol_validate_nontilt.py
+++ b/xmipp3/protocols/protocol_validate_nontilt.py
@@ -341,9 +341,6 @@ def _defineMetadataRootName(self, mdrootname, volId):
             VolPrefix = 'vol%03d_' % (volId)
             return self._getExtraPath(VolPrefix + 'validation.xmd')
 
-    def _definePName(self):
-        fscFn = self._defineMetadataRootName('P')
-        return fscFn
 
     def _defineVolumeName(self, volId):
         fscFn = self._defineMetadataRootName('Volume', volId)
diff --git a/xmipp3/protocols/protocol_volume_local_adjust.py b/xmipp3/protocols/protocol_volume_local_adjust.py
index b166f627b..339f6c1dc 100644
--- a/xmipp3/protocols/protocol_volume_local_adjust.py
+++ b/xmipp3/protocols/protocol_volume_local_adjust.py
@@ -31,12 +31,16 @@
 from pwem.convert import headers
 from pwem.objects import Volume, Transform
 from pwem.protocols import EMProtocol
+from pyworkflow import BETA, UPDATED, NEW, PROD
+
 
 class XmippProtLocalVolAdj(EMProtocol):
     """Protocol to adjust locally volume intensity to a reference volume. Occupancy volume is saved in protocol folder. Based on 
     https://www.sciencedirect.com/science/article/pii/S1047847723000874?via%3Dihub"""
     _label = 'volume local adjustment'
     _possibleOutputs = Volume
+    _devStatus = NEW
+
 
     # --------------------------- DEFINE param functions --------------------------------------------
     def _defineParams(self, form):
diff --git a/xmipp3/viewers/__init__.py b/xmipp3/viewers/__init__.py
index a8dbb1fb5..98ef01723 100644
--- a/xmipp3/viewers/__init__.py
+++ b/xmipp3/viewers/__init__.py
@@ -33,6 +33,7 @@
 from .viewer_deep_consensus import XmippDeepConsensusViewer
 from .viewer_deep_micrograph_cleaner import XmippDeepMicrographViewer
 from .viewer_ml2d import XmippML2DViewer
+from .viewer_movie_alignment import XmippMovieAlignViewer, XmippMovieMaxShiftViewer, XmippMovieAlignConsensusViewer
 #from .viewer_mltomo import XmippMLTomoViewer
 from .viewer_movie_alignment import XmippMovieAlignViewer, XmippMovieMaxShiftViewer
 from .viewer_normalize_strain import XmippNormalizeStrainViewer
@@ -44,11 +45,9 @@
 from .viewer_resolution_deepres import XmippResDeepResViewer
 from .viewer_validate_nontilt import XmippValidateNonTiltViewer
 from .viewer_validate_fscq import XmippProtValFitViewer
-from .viewer_split_volume import XmippViewerSplitVolume
 from .viewer_validate_overfitting import XmippValidateOverfittingViewer
 from .viewer_volume_strain import XmippVolumeStrainViewer
 from .viewer_reconstruct_highres import XmippReconstructHighResViewer
-from .viewer_solid_angles import SolidAnglesViewer
 from .viewer_consensus_classes import XmippConsensusClassesViewer
 from .viewer_extract_asymmetric_unit import viewerXmippProtExtractUnit
 from .viewer_deepEMHancer import viewerXmippProtDeepVolPostProc
@@ -64,7 +63,6 @@
 from .viewer_volume_deform_zernike3d import XmippVolumeDeformZernike3DViewer
 from .viewer_structure_map import XmippProtStructureMapViewer
 from .viewer_pdb_deform_zernike3d import XmippPDBDeformViewer
-from .viewer_metaprotocol_golden_highres import XmippMetaprotocolGoldenHighResViewer
 
 # TODO: Import from continuousflex the modules needed to create clusters (soft dependency)
 from .viewer_analyze_local_ctf import XmippAnalyzeLocalCTFViewer
diff --git a/xmipp3/viewers/viewer.py b/xmipp3/viewers/viewer.py
index f6ede2ae3..961cd7ca7 100644
--- a/xmipp3/viewers/viewer.py
+++ b/xmipp3/viewers/viewer.py
@@ -49,7 +49,6 @@
 from xmipp3.protocols import XmippProtAngularGraphConsistency
 from xmipp3.protocols import XmippProtAssignmentTiltPair
 from xmipp3.protocols import XmippProtMovieGain
-from xmipp3.protocols import XmippProtDeepDenoising
 from .plotter import XmippPlotter
 from xmipp3.protocols import XmippProtTiltAnalysis
 from pwem.viewers.viewers_data import MicrographsView
@@ -76,7 +75,6 @@ class XmippViewer(DataViewer):
                 XmippProtMultiRefAlignability,
                 XmippProtAngularGraphConsistency,
                 XmippProtMovieGain,
-                XmippProtDeepDenoising,
                 XmippProtTiltAnalysis
                 ]
 
@@ -151,26 +149,6 @@ def _visualize(self, obj, **kwargs):
                     self._views.append(xplotter)
 
 
-        elif issubclass(cls, XmippProtDeepDenoising):
-            fn = obj.outputParticles.getFileName()
-            self._views.append(ObjectView(self._project, obj.outputParticles.strId(), fn,
-                                          viewParams={VISIBLE: 'enabled id _filename '
-                                                               '_xmipp_corrDenoisedProjection '
-                                                               '_xmipp_corrDenoisedNoisy '
-                                                               '_xmipp_imageOriginal _xmipp_imageRef',
-                                                      RENDER: '_filename _xmipp_imageOriginal '
-                                                              '_xmipp_imageRef',
-                                                      SORT_BY: 'id',
-                                                      MODE: MODE_MD}))
-
-            md = emlib.MetaData(obj._getExtraPath('particlesDenoised.xmd'))
-            if (md.containsLabel(emlib.MDL_CORR_DENOISED_PROJECTION) and
-                md.containsLabel(emlib.MDL_CORR_DENOISED_NOISY)):
-                xplotter = XmippPlotter(windowTitle="denoised vs proj & denoised vs original")
-                xplotter.createSubPlot("Correlations", "corr_denoised_projection", "corr_denoised_original")
-                xplotter.plotScatterMd(md, mdLabelX=emlib.MDL_CORR_DENOISED_PROJECTION,
-                                mdLabelY=emlib.MDL_CORR_DENOISED_NOISY )
-                self._views.append(xplotter)
 
         elif issubclass(cls, XmippProtMovieGain):
             if obj.hasAttribute('outputGains'):
@@ -225,19 +203,38 @@ def _visualize(self, obj, **kwargs):
 
 
         elif issubclass(cls, XmippProtCompareReprojections):
-                fn = obj.reprojections.getFileName()
-                labels = ('id enabled _index _xmipp_image._filename _xmipp_imageRef._filename '
-                          '_xmipp_imageResidual._filename _xmipp_imageCovariance._filename '
-                          '_xmipp_cost _xmipp_zScoreResCov _xmipp_zScoreResMean _xmipp_zScoreResVar '
-                          '_xmipp_continuousA _xmipp_continuousB _xmipp_continuousX _xmipp_continuousY')
-                labelRender = ("_xmipp_image._filename _xmipp_imageRef._filename "
-                              "_xmipp_imageResidual._filename _xmipp_imageCovariance._filename")
-                self._views.append(ObjectView(self._project, obj.reprojections.strId(), fn,
-                                              viewParams={ORDER: labels,
-                                                          VISIBLE: labels,
-                                                          SORT_BY: '_xmipp_cost asc',
-                                                          RENDER: labelRender,
-                                                          MODE: MODE_MD}))
+            labels = ('id enabled _index _xmipp_image._filename _xmipp_imageRef._filename '
+                      '_xmipp_imageResidual._filename _xmipp_imageCovariance._filename '
+                      '_xmipp_cost _xmipp_zScoreResCov _xmipp_zScoreResMean _xmipp_zScoreResVar '
+                      '_xmipp_continuousA _xmipp_continuousB _xmipp_continuousX _xmipp_continuousY')
+            labelRender = ("_xmipp_image._filename _xmipp_imageRef._filename "
+                           "_xmipp_imageResidual._filename _xmipp_imageCovariance._filename")
+            try:
+                for outputName in obj.readOutputDict().values():
+                    reprojections = getattr(obj, outputName)
+                    fn = reprojections.getFileName()
+                    self._views.append(ObjectView(self._project, reprojections.strId(), fn,
+                                                  viewParams={ORDER: labels,
+                                                              VISIBLE: labels,
+                                                              SORT_BY: '_xmipp_cost asc',
+                                                              RENDER: labelRender,
+                                                              MODE: MODE_MD}))
+            except FileNotFoundError:
+                self._views.append(self.infoMessage('This may be an old version of this protocol, trying older viewer.'))
+                try:
+                    fn = obj.reprojections.getFileName()
+                    self._views.append(ObjectView(self._project, obj.reprojections.strId(), fn,
+                                                  viewParams={ORDER: labels,
+                                                              VISIBLE: labels,
+                                                              SORT_BY: '_xmipp_cost asc',
+                                                              RENDER: labelRender,
+                                                              MODE: MODE_MD}))
+                except Exception as e:
+                    self._views.append(self.infoMessage("There is a problem with the output: %s" % e))
+
+            except Exception as e:
+                self._views.append(self.infoMessage("There is not an output yet: %s" % e))
+
 
         elif issubclass(cls, XmippProtCompareAngles):
                 fn = obj.outputParticles.getFileName()
diff --git a/xmipp3/viewers/viewer_extract_asymmetric_unit.py b/xmipp3/viewers/viewer_extract_asymmetric_unit.py
index 08d58c0a9..1f4e615d2 100644
--- a/xmipp3/viewers/viewer_extract_asymmetric_unit.py
+++ b/xmipp3/viewers/viewer_extract_asymmetric_unit.py
@@ -141,8 +141,6 @@ def _showVolumesChimera(self):
         if innerRadius < 0:
            innerRadius = 0
         d['innerRadius'] = innerRadius * sampling
-
-        d['innerRadius'] = self.protocol.innerRadius.get() * sampling
         d['symmetry'] = Chimera.getSymmetry(XMIPP_TO_SCIPION[self.protocol.symmetryGroup.get()])
 
         if self.protocol.symmetryGroup >= XMIPP_I222:
diff --git a/xmipp3/viewers/viewer_movie_alignment.py b/xmipp3/viewers/viewer_movie_alignment.py
index d40a8e504..80234caf4 100644
--- a/xmipp3/viewers/viewer_movie_alignment.py
+++ b/xmipp3/viewers/viewer_movie_alignment.py
@@ -29,6 +29,7 @@
 from pwem.viewers import showj, EmProtocolViewer, ObjectView
 from pwem.objects import SetOfMicrographs, SetOfMovies
 
+from xmipp3.protocols.protocol_movie_alignment_consensus import XmippProtConsensusMovieAlignment
 from xmipp3.protocols.protocol_movie_opticalflow import (XmippProtOFAlignment,
                                                  OBJCMD_MOVIE_ALIGNCARTESIAN)
 from xmipp3.protocols.protocol_flexalign import XmippProtFlexAlign
@@ -107,6 +108,7 @@ def _visualizeMics(self, e=None):
     def _visualizeMicsDiscarded(self, e=None):
         return self._visualizeAny(lambda x: x.endswith('Discarded'))
 
+
 def getViewParams():
     plotLabels = ('psdCorr._filename plotPolar._filename '
                   'plotCart._filename plotGlobal._filename')
@@ -120,3 +122,44 @@ def getViewParams():
                   }
 
     return viewParams
+
+
+class XmippMovieAlignConsensusViewer(Viewer):
+    _targets = [XmippProtConsensusMovieAlignment]
+    _environments = [DESKTOP_TKINTER, WEB_DJANGO]
+
+    _label = 'viewer Global Alignment Consensus'
+
+    def _visualize(self, obj, **kwargs):
+        views = []
+
+        if obj.hasAttribute('outputMicrographsDiscarded'):
+            views.append(ObjectView(self._project, obj.strId(),
+                                    obj.outputMicrographsDiscarded.getFileName(),
+                                    viewParams=getViewParamsConsensus()))
+        if obj.hasAttribute('outputMicrographs'):
+            views.append(ObjectView(self._project, obj.strId(),
+                                    obj.outputMicrographs.getFileName(),
+                                    viewParams=getViewParamsConsensus()))
+
+        if not(obj.hasAttribute('outputMicrographs')) and not(obj.hasAttribute('outputMicrographsDiscarded')):
+            views.append(self.infoMessage("Output (micrographs or movies) has "
+                                          "not been produced yet."))
+
+        return views
+
+
+def getViewParamsConsensus():
+    plotLabels = ('psdCorr._filename plotPolar._filename '
+                  'plotCart._filename plotGlobal._filename')
+
+    labels = plotLabels + ' _filename ' + '_alignment_corr ' + '_alignment_rmse_error '
+    viewParams = {showj.MODE: showj.MODE_MD,
+                  showj.ORDER: 'id ' + labels,
+                  showj.VISIBLE: 'id ' + labels,
+                  showj.RENDER: plotLabels,
+                  showj.ZOOM: 10,
+                  showj.OBJCMDS: "'%s'" % OBJCMD_MOVIE_ALIGNCARTESIAN
+                  }
+
+    return viewParams
diff --git a/xmipp3/viewers/viewer_projmatch.py b/xmipp3/viewers/viewer_projmatch.py
index 57caee6f4..6dc546f9b 100644
--- a/xmipp3/viewers/viewer_projmatch.py
+++ b/xmipp3/viewers/viewer_projmatch.py
@@ -783,7 +783,7 @@ def _getGuinierValue(self, guinierFn, col):
 #===============================================================================
 # Utils Functions
 #===============================================================================
-    def createScipionPartView(self, partSet, viewParams={}):
+    def createScipionPartView(self, partSet):
         inputParticlesId = self.protocol.inputParticles.get().strId()
         filename = partSet.getFileName()
         
diff --git a/xmipp3/viewers/viewer_structure_map.py b/xmipp3/viewers/viewer_structure_map.py
index a1af5054d..1fb81b6d8 100644
--- a/xmipp3/viewers/viewer_structure_map.py
+++ b/xmipp3/viewers/viewer_structure_map.py
@@ -28,7 +28,7 @@
 from scipy import ndimage
 from scipy.spatial import KDTree
 import tkinter as tk
-matplotlib.use("TkAgg")
+#matplotlib.use("TkAgg")
 import matplotlib.pyplot as plt
 from mpl_toolkits.mplot3d import Axes3D
 from matplotlib.widgets import RadioButtons
diff --git a/xmipp3/viewers/viewer_validate_fscq.py b/xmipp3/viewers/viewer_validate_fscq.py
index 9ca459510..0ffce413e 100644
--- a/xmipp3/viewers/viewer_validate_fscq.py
+++ b/xmipp3/viewers/viewer_validate_fscq.py
@@ -29,23 +29,15 @@
 from tkinter import Tk, CENTER, Scrollbar, TOP, BOTH, Y
 
 import numpy as np
-import matplotlib.pyplot as plt
-from matplotlib import cm
-import matplotlib.colors as mcolors
-from pyworkflow.utils import getExt, removeExt
-from os.path import abspath
 from pwem.viewers import (LocalResolutionViewer, EmPlotter, ChimeraView,
-                          DataView, ChimeraAttributeViewer)
+                          ChimeraAttributeViewer)
 from pyworkflow.gui import *
 from pyworkflow.protocol.params import (LabelParam, EnumParam)
 from pyworkflow.viewer import ProtocolViewer, DESKTOP_TKINTER
 from xmipp3.protocols.protocol_validate_fscq import (XmippProtValFit, 
                                                      RESTA_FILE_MRC, 
                                                      OUTPUT_PDBMRC_FILE, 
-                                                     PDB_VALUE_FILE,
-                                                     RESTA_FILE_NORM,
-                                                     PDB_NORM_FILE)
-
+                                                     RESTA_FILE_NORM)
 
 class XmippProtValFitViewer(LocalResolutionViewer, ChimeraAttributeViewer):
     """
@@ -77,7 +69,6 @@ def _defineParams(self, form):
                important=True,
                label='Display FSC-Qr Volume Output')
 
-
         group.addParam('displayPDB', EnumParam,
                       choices=['by residue', 'by atom'],
                       default=0, important=True,
@@ -132,10 +123,14 @@ def _create_legend(self, scale):
         f = open(fnCmd, 'w')
         f.write("from chimerax.core.commands import run\n")
         f.write("from chimerax.graphics.windowsize import window_size\n")
-        f.write("from PyQt5.QtGui import QFontMetrics\n")
-        f.write("from PyQt5.QtGui import QFont\n")
+        f.write("try:\n")
+        f.write("    from PyQt5.QtGui import QFontMetrics\n")
+        f.write("    from PyQt5.QtGui import QFont\n")
+        f.write("except ModuleNotFoundError:\n")
+        f.write("    from PyQt6.QtGui import QFontMetrics\n")
+        f.write("    from PyQt6.QtGui import QFont\n")
+
         f.write("run(session, 'set bgColor white')\n")
-                
         # get window size so we can place labels properly
         f.write("v = session.main_view\n")
         f.write("vx,vy=v.window_size\n")
@@ -160,7 +155,7 @@ def _create_legend(self, scale):
 
     def _visualize_vol(self, obj, **args):
         
-        self._create_legend(3);
+        self._create_legend(3)
         
         fnCmd = self.protocol._getExtraPath("chimera_output.py")
         f = open(fnCmd, 'a')
@@ -183,7 +178,7 @@ def _visualize_vol(self, obj, **args):
     
     def _visualize_norm_vol(self, obj, **args):
         
-        self._create_legend(1.5);
+        self._create_legend(1.5)
         
         fnCmd = self.protocol._getExtraPath("chimera_output.py")
         f = open(fnCmd, 'a')
@@ -206,12 +201,12 @@ def _visualize_norm_vol(self, obj, **args):
     
     def _visualize_pdb(self, obj, **args):
         
-        self._create_legend(3);
+        self._create_legend(3)
         
         fnCmd = self.protocol._getExtraPath("chimera_output.py")
         f = open(fnCmd, 'a')
                
-        f.write("run(session, 'open %s')\n" % self.protocol._getFileName(PDB_VALUE_FILE))
+        f.write("run(session, 'open %s')\n" % self.protocol.getFSCQFile())
          
         if self.displayPDB == self.RESIDUE:
             f.write("run(session, 'cartoon')\n")
@@ -234,12 +229,12 @@ def _visualize_pdb(self, obj, **args):
     
     def _visualize_norm_pdb(self, obj, **args):
         
-        self._create_legend(1.5);
+        self._create_legend(1.5)
         
         fnCmd = self.protocol._getExtraPath("chimera_output.py")
         f = open(fnCmd, 'a')
                
-        f.write("run(session, 'open %s')\n" % self.protocol._getFileName(PDB_NORM_FILE))
+        f.write("run(session, 'open %s')\n" % self.protocol.getNormFSCQFile())
          
         if self.displayNormPDB == self.RESIDUE:
             f.write("run(session, 'cartoon')\n")
@@ -262,62 +257,56 @@ def _visualize_norm_pdb(self, obj, **args):
  
     def _calculate_fscq(self, obj, **args):
 
-        fnRoot = os.path.abspath(self.protocol._getExtraPath())
-        bool = 0
-        overfitting_list = []
-        poorfitting_list = []
-        with open(fnRoot + '/' + PDB_VALUE_FILE) as f:
-
-            lines_data = f.readlines()
-            for j, lin in enumerate(lines_data):
+        status = False
+        overfittingList = []
+        poorfittingList = []
+        with open(os.path.abspath(self.protocol.getFSCQFile())) as f:
+            linesData = f.readlines()
+            for j, lin in enumerate(linesData):
                 if lin.startswith('ATOM') or lin.startswith('HETATM'):
                     resnumber = int(lin[22:26])
 
-                    if bool == 1 and resnumber == resnumber_ctl:
-                        resatomname = lin[12:16].strip()
+                    if status and resnumber == resnumberCtl:
                         resname = lin[17:20].strip()
                         chain = lin[21]
                         fscq = float(lin[54:60])
-                        current_frag.append(fscq)
+                        currentFrag.append(fscq)
 
-                    elif bool == 1 and resnumber != resnumber_ctl:
-                        meanFscq = np.mean(current_frag)
-                        current_frag.sort()
+                    elif status and resnumber != resnumberCtl:
+                        meanFscq = np.mean(currentFrag)
+                        currentFrag.sort()
 
-                        if current_frag[0] <= -1:
-                            overfitting_list.append((resname, resnumber_ctl, chain,
-                                                     current_frag[0], meanFscq))
+                        if currentFrag[0] <= -1:
+                            overfittingList.append((resname, resnumberCtl, chain,
+                                                     currentFrag[0], meanFscq))
 
-                        if current_frag[-1] >= 1:
-                            poorfitting_list.append((resname, resnumber_ctl, chain,
-                                                     current_frag[-1],
+                        if currentFrag[-1] >= 1:
+                            poorfittingList.append((resname, resnumberCtl, chain,
+                                                     currentFrag[-1],
                                                      meanFscq))
                      
-                        current_frag = []
-                        resnumber_ctl = resnumber
-                        resatomname = lin[12:16].strip()
+                        currentFrag = []
+                        resnumberCtl = resnumber
                         resname = lin[17:20].strip()
                         chain = lin[21]
                         fscq = float(lin[54:60])
 
-                        current_frag.append(fscq)
+                        currentFrag.append(fscq)
 
                     else:
-                        bool = 1
-                        current_frag = []
-                        lines_aa = []
-                        resnumber_ctl = resnumber
-                        resatomname = lin[12:16].strip()
+                        status = True
+                        currentFrag = []
+                        resnumberCtl = resnumber
                         resname = lin[17:20].strip()
                         chain = lin[21]
                         fscq = float(lin[54:60])
 
-                        current_frag.append(fscq)
+                        currentFrag.append(fscq)
                         
         if obj == 'calculateFscqNeg':
-            return overfitting_list
+            return overfittingList
         else:
-            return poorfitting_list
+            return poorfittingList
 
     def _statistics(self, obj, **args):
         mainFrame = Tk()


Test patch:
diff --git a/xmipp3/tests/test_protocol_deep_denoising.py b/xmipp3/legacy/tests/test_protocol_deep_denoising.py
similarity index 100%
rename from xmipp3/tests/test_protocol_deep_denoising.py
rename to xmipp3/legacy/tests/test_protocol_deep_denoising.py
diff --git a/xmipp3/tests/test_protocols_metaprotocol_golden_highres.py b/xmipp3/legacy/tests/test_protocols_metaprotocol_golden_highres.py
similarity index 100%
rename from xmipp3/tests/test_protocols_metaprotocol_golden_highres.py
rename to xmipp3/legacy/tests/test_protocols_metaprotocol_golden_highres.py
diff --git a/xmipp3/tests/test_protocols_solid_angles.py b/xmipp3/legacy/tests/test_protocols_solid_angles.py
similarity index 100%
rename from xmipp3/tests/test_protocols_solid_angles.py
rename to xmipp3/legacy/tests/test_protocols_solid_angles.py
diff --git a/xmipp3/tests/test_convert_xmipp.py b/xmipp3/tests/test_convert_xmipp.py
index 68dce66c6..56c4f48f3 100644
--- a/xmipp3/tests/test_convert_xmipp.py
+++ b/xmipp3/tests/test_convert_xmipp.py
@@ -1017,45 +1017,45 @@ def test_writeSetOfDefocusGroups(self):
         #TODO: FIX THIS test according to the new SetOfDefocusGroup
         return
         #reference metadata
-        md = emlib.MetaData()
-        objId = md.addObject()
-        defocusGroupRow = metadata.Row()
-
-        defocusGroupRow.setValue(emlib.MDL_ENABLED, 1)
-        defocusGroupRow.setValue(emlib.MDL_CTF_GROUP, 1)
-        defocusGroupRow.setValue(emlib.MDL_MIN, 2000.)
-        defocusGroupRow.setValue(emlib.MDL_MAX, 2500.)
-        defocusGroupRow.setValue(emlib.MDL_AVG, 2100.)
-        defocusGroupRow.writeToMd(md, objId)
-
-        objId = md.addObject()
-        defocusGroupRow.setValue(emlib.MDL_ENABLED, 1)
-        defocusGroupRow.setValue(emlib.MDL_CTF_GROUP, 2)
-        defocusGroupRow.setValue(emlib.MDL_MIN, 3000.)
-        defocusGroupRow.setValue(emlib.MDL_MAX, 5500.)
-        defocusGroupRow.setValue(emlib.MDL_AVG, 5000.)
-        defocusGroupRow.writeToMd(md, objId)
+        # md = emlib.MetaData()
+        # objId = md.addObject()
+        # defocusGroupRow = metadata.Row()
+        #
+        # defocusGroupRow.setValue(emlib.MDL_ENABLED, 1)
+        # defocusGroupRow.setValue(emlib.MDL_CTF_GROUP, 1)
+        # defocusGroupRow.setValue(emlib.MDL_MIN, 2000.)
+        # defocusGroupRow.setValue(emlib.MDL_MAX, 2500.)
+        # defocusGroupRow.setValue(emlib.MDL_AVG, 2100.)
+        # defocusGroupRow.writeToMd(md, objId)
+        #
+        # objId = md.addObject()
+        # defocusGroupRow.setValue(emlib.MDL_ENABLED, 1)
+        # defocusGroupRow.setValue(emlib.MDL_CTF_GROUP, 2)
+        # defocusGroupRow.setValue(emlib.MDL_MIN, 3000.)
+        # defocusGroupRow.setValue(emlib.MDL_MAX, 5500.)
+        # defocusGroupRow.setValue(emlib.MDL_AVG, 5000.)
+        # defocusGroupRow.writeToMd(md, objId)
+        # #
+        # fnScipion = self.getOutputPath("writeSetOfDefocusGroups.sqlite")
+        # setOfDefocus = SetOfDefocusGroup(filename=fnScipion)
+        #
+        # df = DefocusGroup()
+        # df.setDefocusMin(2000.)
+        # df.setDefocusMax(2500.)
+        # df.setDefocusAvg(2100.)
+        # setOfDefocus.append(df)
+        #
+        # df.cleanObjId()
+        # df.setDefocusMin(3000)
+        # df.setDefocusMax(5500)
+        # df.setDefocusAvg(5000)
+        # setOfDefocus.append(df)
+        #
+        # fnXmipp = self.getOutputPath("writeSetOfDefocusGroups.xmd")
+        # fnScipion = self.getOutputPath("writeSetOfDefocusGroups2.xmd")
+        # writeSetOfDefocusGroups(setOfDefocus, fnXmipp)
+        # mdAux = emlib.MetaData(fnXmipp)
+        # md.write(fnScipion)
+        # print("Comparing metadatas: \n%s\n%s" % (fnXmipp, fnScipion))
+        # self.assertEqual(md, mdAux, "test writeSetOfDefocusGroups fails")
         #
-        fnScipion = self.getOutputPath("writeSetOfDefocusGroups.sqlite")
-        setOfDefocus = SetOfDefocusGroup(filename=fnScipion)
-
-        df = DefocusGroup()
-        df.setDefocusMin(2000.)
-        df.setDefocusMax(2500.)
-        df.setDefocusAvg(2100.)
-        setOfDefocus.append(df)
-        
-        df.cleanObjId()
-        df.setDefocusMin(3000)
-        df.setDefocusMax(5500)
-        df.setDefocusAvg(5000)
-        setOfDefocus.append(df)
-        
-        fnXmipp = self.getOutputPath("writeSetOfDefocusGroups.xmd")
-        fnScipion = self.getOutputPath("writeSetOfDefocusGroups2.xmd")
-        writeSetOfDefocusGroups(setOfDefocus, fnXmipp)
-        mdAux = emlib.MetaData(fnXmipp)
-        md.write(fnScipion)
-        print("Comparing metadatas: \n%s\n%s" % (fnXmipp, fnScipion))
-        self.assertEqual(md, mdAux, "test writeSetOfDefocusGroups fails")
-       
diff --git a/xmipp3/tests/test_protocol_deep_center_assignment.py b/xmipp3/tests/test_protocol_deep_center_assignment.py
deleted file mode 100644
index 700b1d7b2..000000000
--- a/xmipp3/tests/test_protocol_deep_center_assignment.py
+++ /dev/null
@@ -1,132 +0,0 @@
-# **************************************************************************
-# *
-# * Authors:    Adrian Sansinena Rodriguez  (adrian.sansinena@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'
-# *
-# **************************************************************************
-
-################## QUALITY DISCLAIMER ##################
-# This tests only proves that this protocol can successfully run without errors!
-# It does not guarantee at all that the results will be any good.
-# The input train set is the same as the set to predict, which causes artificially great results,
-# that will get worse once real use cases are introduced. However, changing that makes no sense 
-# either, since the quality of the results cannot be checked, due to the non-deterministic nature
-# of Neural Networks, specially regarding their training.
-
-# Scipion em imports
-from pwem.protocols import ProtImportParticles, ProtSubSet, exists
-from pyworkflow.tests import BaseTest, DataSet, setupTestProject
-
-# Plugin imports
-from ..protocols import XmippProtDeepGlobalAssignment, XmippProtDeepCenter
-
-# ---------------------------- COMMON FUNCTIONS & VARIABLES ----------------------------
-genericError = "There was a problem while running this protocol."
-
-def runImportParticles(cls):
-    """ Import Particles. """
-    args = {'importFrom': ProtImportParticles.IMPORT_FROM_SCIPION,
-            'sqliteFile': cls.particles,
-            'amplitudConstrast': 0.1,
-            'sphericalAberration': 2.0,
-            'voltage': 200,
-            'samplingRate': 0.99,
-            'haveDataBeenPhaseFlipped': True
-            }
-
-    protImport = cls.newProtocol(ProtImportParticles, **args)
-    protImport.setObjLabel('import particles')
-    cls.launchProtocol(protImport)
-    return protImport
-
-def setUpClassBase(cls):
-    setupTestProject(cls)
-
-    # Data
-    cls.dataset = DataSet.getDataSet('10010')
-    cls.particles = cls.dataset.getFile('particles')
-    cls.protImportParts = runImportParticles(cls)
-
-class TestDeepGlobalAssignment(BaseTest):
-    @classmethod
-    def setUpClass(cls):
-        setUpClassBase(cls)
-
-    def test(self):
-        # Creating subset of particles (same for train and test)
-        subset = self.newProtocol(ProtSubSet,
-                                  inputFullSet=self.protImportParts.outputParticles,
-                                  chooseAtRandom=True,
-                                  nElements=400)
-        self.launchProtocol(subset)
-
-        deepGA = self.newProtocol(XmippProtDeepGlobalAssignment,
-                                  inputImageSet=subset.outputParticles,
-                                  inputTrainSet=subset.outputParticles,
-                                  Xdim=128,
-                                  numModels=5,
-                                  numEpochs=1,
-                                  batchSize=32,
-                                  learningRate=0.001,
-                                  sigma=8,
-                                  patience=5)
-        self.launchProtocol(deepGA)
-        self.assertIsNotNone(deepGA.outputParticles, genericError)
-
-class TestDeepCenter(BaseTest):
-    @classmethod
-    def setUpClass(cls):
-        setUpClassBase(cls)
-
-    def test(self):
-        subset = self.newProtocol(ProtSubSet,
-                                  inputFullSet=self.protImportParts.outputParticles,
-                                  chooseAtRandom=True,
-                                  nElements=400)
-        self.launchProtocol(subset)
-
-        deepCenter = self.newProtocol(XmippProtDeepCenter,
-                                   inputImageSet=subset.outputParticles,
-                                   trainModels=True,
-                                   inputTrainSet=subset.outputParticles,
-                                   numModels=5,
-                                   numEpochs=1,
-                                   batchSize=32,
-                                   learningRate=0.001,
-                                   sigma=8,
-                                   patience=5)
-        self.launchProtocol(deepCenter)
-
-        fnModel0 = deepCenter._getExtraPath("model0.h5")
-        fnModel2 = deepCenter._getExtraPath("model2.h5")
-
-        self.assertTrue(exists(fnModel0), fnModel0 + " does not exist")
-        self.assertTrue(exists(fnModel2), fnModel2 + " does not exist")
-
-        self.assertIsNotNone(deepCenter.outputParticles, genericError)
-
-        deepCenter2 = self.newProtocol(XmippProtDeepCenter,
-                                       inputImageSet=subset.outputParticles,
-                                       trainModels=False)
-        self.launchProtocol(deepCenter2)
-
-        self.assertIsNotNone(deepCenter2.outputParticles, genericError)
diff --git a/xmipp3/tests/test_protocols_xmipp_2d.py b/xmipp3/tests/test_protocols_xmipp_2d.py
index 9d4241d33..1168bebb4 100644
--- a/xmipp3/tests/test_protocols_xmipp_2d.py
+++ b/xmipp3/tests/test_protocols_xmipp_2d.py
@@ -28,6 +28,8 @@
 
 from __future__ import print_function
 
+import os.path
+
 from pyworkflow.tests.test_utils import wait
 from pyworkflow.utils import greenStr, magentaStr
 from pyworkflow.plugin import Domain
@@ -917,7 +919,8 @@ def setUpClass(cls):
         setupTestProject(cls)
         TestXmippBase.setData('mda')
         cls.protImport = cls.runImportParticles(cls.particlesFn, 3.5)
-        cls.protImportAvgs = cls.runImportAverages(cls.particlesDir + '/img00007[1-4].spi', 3.5)
+        cls.protImportAvgs = cls.runImportAverages(os.path.dirname(cls.particlesFn) +
+                                                   '/img00007[1-4].spi', 3.5)
     
     def test_cl2d(self):
         print("Run CL2D")
@@ -1174,9 +1177,12 @@ class TestXmippCompareReprojections(TestXmippBase):
     def setUpClass(cls):
         setupTestProject(cls)
         TestXmippBase.setData('mda')
+        cls.dsRelion = DataSet.getDataSet('relion_tutorial')
         cls.protImportPart = cls.runImportParticles(cls.particlesFn, 3.5)
         cls.protImportAvgs = cls.runImportAverages(cls.particlesFn, 3.5)
         cls.protImportVol = cls.runImportVolume(cls.volumesFn, 3.5)
+        cls.protImport2D = cls.importFromRelion2D()
+        cls.protImport3D = cls.importFromRelionRefine3D()
         cls.protCreateMask = cls.newProtocol(XmippProtCreateMask3D, threshold=0.01)
         cls.protCreateMask.inputVolume.set(cls.protImportVol.outputVolume)
         cls.launchProtocol(cls.protCreateMask)
@@ -1193,74 +1199,120 @@ def setUpClass(cls):
         cls.protProjMatch.inputParticles.set(cls.protImportAvgs.outputAverages)
         cls.protProjMatch.input3DReferences.set(cls.protImportVol.outputVolume)
         cls.launchProtocol(cls.protProjMatch)
-    
+
+    @classmethod
+    def importFromRelionRefine3D(cls):
+        """ Import aligned Particles
+        """
+        protImport3D = cls.newProtocol(emprot.ProtImportParticles,
+                                           objLabel='particles from relion (auto-refine 3d)',
+                                           importFrom=emprot.ProtImportParticles.IMPORT_FROM_RELION,
+                                           starFile=
+                                           cls.dsRelion.getFile('import/classify3d/extra/relion_it015_data.star'),
+                                           magnification=10000,
+                                           samplingRate=7.08,
+                                           haveDataBeenPhaseFlipped=True)
+        cls.launchProtocol(protImport3D)
+        return protImport3D
+
+    @classmethod
+    def importFromRelion2D(cls):
+        """
+        Import an EMX file with Particles and defocus
+        """
+        protImport2D = cls.newProtocol(emprot.ProtImportParticles,
+                                           objLabel='from relion (classify 2d)',
+                                           importFrom=emprot.ProtImportParticles.IMPORT_FROM_RELION,
+                                           starFile=cls.dsRelion.getFile(
+                                             'import/classify2d/extra/relion_it015_data.star'),
+                                           magnification=10000,
+                                           samplingRate=7.08,
+                                           haveDataBeenPhaseFlipped=True)
+        cls.launchProtocol(protImport2D)
+        return protImport2D
+
     def test_particles1(self):
         print("Run Compare Reprojections from classes")
         prot = self.newProtocol(XmippProtCompareReprojections,
-                                        symmetryGroup="d6", numberOfMpi=5)
-        prot.inputSet.set(self.protClassify.outputClasses)
-        prot.inputVolume.set(self.protImportVol.outputVolume)
+                                        symmetryGroup="d6", numberOfMpi=5, doRanking=False)
+        prot.inputSet2D.set(self.protClassify.outputClasses)
+        prot.inputSet3D.set(self.protImportVol.outputVolume)
         self.launchProtocol(prot)
-        self.assertIsNotNone(prot.reprojections, "There was a problem with Compare Reprojections from classes")
+        self.assertIsNotNone(prot.reprojections_vol200, "There was a problem with Compare Reprojections from classes")
 
     def test_particles2(self):
         print("Run Compare Reprojections from averages")
         prot = self.newProtocol(XmippProtCompareReprojections,
-                                        symmetryGroup="d6", numberOfMpi=5)
-        prot.inputSet.set(self.protImportAvgs.outputAverages)
-        prot.inputVolume.set(self.protImportVol.outputVolume)
+                                        symmetryGroup="d6", numberOfMpi=5, doRanking=False)
+        prot.inputSet2D.set(self.protImportAvgs.outputAverages)
+        prot.inputSet3D.set(self.protImportVol.outputVolume)
         self.launchProtocol(prot)
-        self.assertIsNotNone(prot.reprojections, "There was a problem with Compare Reprojections from averages")
+        self.assertIsNotNone(prot.reprojections_vol200, "There was a problem with Compare Reprojections from averages")
 
     def test_particles3(self):
         print("Run Compare Reprojections from projections with angles")
         prot = self.newProtocol(XmippProtCompareReprojections,
-                                        symmetryGroup="d6", numberOfMpi=5)
-        prot.inputSet.set(self.protProjMatch.outputParticles)
-        prot.inputVolume.set(self.protImportVol.outputVolume)
+                                        symmetryGroup="d6", numberOfMpi=5, doRanking=False)
+        prot.inputSet2D.set(self.protProjMatch.outputParticles)
+        prot.inputSet3D.set(self.protImportVol.outputVolume)
         self.launchProtocol(prot)
-        self.assertIsNotNone(prot.reprojections,
+        self.assertIsNotNone(prot.reprojections_vol200,
                              "There was a problem with Compare Reprojections from projections with angles")
 
     def test_particles4(self):
         print("Run Compare Reprojections from classes evaluating residuals")
         prot = self.newProtocol(XmippProtCompareReprojections,
-                                symmetryGroup="d6", numberOfMpi=5, doEvaluateResiduals=True)
-        prot.inputSet.set(self.protClassify.outputClasses)
-        prot.inputVolume.set(self.protImportVol.outputVolume)
+                                symmetryGroup="d6", numberOfMpi=5, doEvaluateResiduals=True, doRanking=False)
+        prot.inputSet2D.set(self.protClassify.outputClasses)
+        prot.inputSet3D.set(self.protImportVol.outputVolume)
         self.launchProtocol(prot)
-        self.assertIsNotNone(prot.reprojections,
+        self.assertIsNotNone(prot.reprojections_vol200,
                              "There was a problem with Compare Reprojections from classes evaluating residuals")
 
     def test_particles5(self):
         print("Run Compare Reprojections from averages evaluating residuals")
         prot = self.newProtocol(XmippProtCompareReprojections,
-                                symmetryGroup="d6", numberOfMpi=5, doEvaluateResiduals=True)
-        prot.inputSet.set(self.protImportAvgs.outputAverages)
-        prot.inputVolume.set(self.protImportVol.outputVolume)
+                                symmetryGroup="d6", numberOfMpi=5, doEvaluateResiduals=True, doRanking=False)
+        prot.inputSet2D.set(self.protImportAvgs.outputAverages)
+        prot.inputSet3D.set(self.protImportVol.outputVolume)
         self.launchProtocol(prot)
-        self.assertIsNotNone(prot.reprojections,
+        self.assertIsNotNone(prot.reprojections_vol200,
                              "There was a problem with Compare Reprojections from averages evaluating residuals")
 
     def test_particles6(self):
         print("Run Compare Reprojections from projections with angles evaluating residuals")
         prot = self.newProtocol(XmippProtCompareReprojections,
-                                symmetryGroup="d6", numberOfMpi=5, doEvaluateResiduals=True)
-        prot.inputSet.set(self.protProjMatch.outputParticles)
-        prot.inputVolume.set(self.protImportVol.outputVolume)
+                                symmetryGroup="d6", numberOfMpi=5, doEvaluateResiduals=True, doRanking=False)
+        prot.inputSet2D.set(self.protProjMatch.outputParticles)
+        prot.inputSet3D.set(self.protImportVol.outputVolume)
         self.launchProtocol(prot)
-        self.assertIsNotNone(prot.reprojections, "There was a problem with Compare Reprojections from projections"
-                                                 " with angles evaluating residuals")
+        self.assertIsNotNone(prot.reprojections_vol200,
+                             "There was a problem with Compare Reprojections from projections"
+                             " with angles evaluating residuals")
 
     def test_particles7(self):
         print("Run Compare Reprojections from projections with angles evaluating residuals without user mask")
         prot = self.newProtocol(XmippProtCompareReprojections,
-                                symmetryGroup="d6", numberOfMpi=5, doEvaluateResiduals=True)
-        prot.inputSet.set(self.protProjMatch.outputParticles)
-        prot.inputVolume.set(self.protImportVol.outputVolume)
+                                symmetryGroup="d6", numberOfMpi=5, doEvaluateResiduals=True, doRanking=False)
+        prot.inputSet2D.set(self.protProjMatch.outputParticles)
+        prot.inputSet3D.set(self.protImportVol.outputVolume)
+        self.launchProtocol(prot)
+        self.assertIsNotNone(prot.reprojections_vol200,
+                             "There was a problem with Compare Reprojections from projections"
+                             " with angles evaluating residuals without user mask")
+
+    def test_ranking3D(self):
+        print("Run Compare Reprojections from 3d and 2d classes and then make a ranking with the best volume")
+        prot = self.newProtocol(XmippProtCompareReprojections,
+                                symmetryGroup="d6", numberOfMpi=5)
+        prot.inputSet2D.set(self.protImport2D.outputClasses)
+        prot.inputSet3D.set(self.protImport3D.outputClasses)
         self.launchProtocol(prot)
-        self.assertIsNotNone(prot.reprojections, "There was a problem with Compare Reprojections from projections"
-                                                 " with angles evaluating residuals without user mask")
+        self.assertIsNotNone(prot.reprojections_vol1, "There was a problem with Compare Reprojections vol 1")
+        self.assertIsNotNone(prot.reprojections_vol2, "There was a problem with Compare Reprojections vol 2")
+        self.assertIsNotNone(prot.reprojections_vol3, "There was a problem with Compare Reprojections vol 3")
+        self.assertSetSize(prot.particles_bestVol, size=899)
+        self.assertIsNotNone(prot.bestVolume, "There is a problem with the ranking the best volume does not exist")
 
 class TestXmippCreateGallery(TestXmippBase):
     """This class check if the protocol create gallery in Xmipp works properly."""
@@ -1450,6 +1502,105 @@ def testXmippScreenDeepLearning(self):
         # Check the sampling rate of the first particle
         self.assertEqual(protScreenDeepLearning.outputParticles.getFirstItem().getSamplingRate(), 7.08, (MSG_WRONG_SAMPLING, "particles"))
 
+class TestXmippClassifyPca(TestXmippBase):
+    """This class check if the protocol Classify PCA (static and in streaming) in Xmipp works properly."""
+
+    @classmethod
+    def setUpClass(cls):
+        setupTestProject(cls)
+        TestXmippBase.setData('mda')
+        cls.dsRelion = DataSet.getDataSet('relion_tutorial')
+        cls.protImport = cls.importParticles()
+
+    @classmethod
+    def importParticles(cls):
+        pathToFile = 'import/case2/relion_it015_data.star'
+        importProt = cls.newProtocol(emprot.ProtImportParticles,
+                                     objLabel='from relion (auto-refine 3d)',
+                                     importFrom=emprot.ProtImportParticles.IMPORT_FROM_RELION,
+                                     starFile=cls.dsRelion.getFile(pathToFile),
+                                     magnification=10000,
+                                     samplingRate=7.08,
+                                     haveDataBeenPhaseFlipped=True
+                                     )
+        cls.launchProtocol(importProt)
+
+        return importProt
+
+    def _updateProtocol(self, prot):
+        prot2 = getProtocolFromDb(prot.getProject().path,
+                                  prot.getDbPath(),
+                                  prot.getObjId())
+        # Close DB connections
+        prot2.getProject().closeMapper()
+        prot2.closeMappers()
+        return prot2
+
+    def test_ClassifyPCAStreaming(self):
+        print("Run 1st Classify PCA Static")
+        protPCA1 = self.newProtocol(XmippProtClassifyPca,
+                                    objLabel="Classify Pca - static version",
+                                    numberOfClasses=5, numberOfMpi=4, numberOfThreads=1)
+        protPCA1.inputParticles.set(self.protImport.outputParticles)
+        self.proj.scheduleProtocol(protPCA1)
+
+        print("Start Streaming Particles")
+        protStream = self.newProtocol(emprot.ProtCreateStreamData, setof=3,
+                                      creationInterval=5, samplingRate=7.08, nDim=6000, groups=500)
+        protStream.inputParticles.set(self.protImport.outputParticles)
+        self.proj.scheduleProtocol(protStream, prerequisites=[protPCA1.getObjId()])
+
+        print("Run 2nd Classify PCA Streaming")
+        protPCA2 = self.newProtocol(XmippProtClassifyPcaStreaming,
+                                    objLabel="Classify Pca streaming - update classes",
+                                    training=2000,
+                                    correctCtf=False,
+                                    mode=XmippProtClassifyPcaStreaming.UPDATE_CLASSES
+                                    )
+        protPCA2.initialClasses.set(protPCA1)
+        protPCA2.initialClasses.setExtended("outputAverages")
+
+        protPCA2.inputParticles.set(protStream)
+        protPCA2.inputParticles.setExtended("outputParticles")
+
+        self.proj.scheduleProtocol(protPCA2)
+
+        # when the first 2D Classification (static mode) finishes must have 5 classes.
+        self.checkResults(protPCA1, 5, "Static mode classification")
+        # when the 2D update Classification (streaming mode) finishes must have 5 classes.
+        self.checkResults(protPCA2, 5, "Streaming mode update initial classification")
+
+        time.sleep(5)  # sometimes is not ready yet
+
+        # Check a final update
+        self.checkFinal2DClasses1st(protPCA1, msg="Initial classification mode fails")
+        self.checkFinal2DClasses2nd(protPCA2, msg="Update classification in streaming mode fails")
+
+    def checkFinal2DClasses1st(self, prot, outputName="outputClasses", msg=''):
+        prot = self._updateProtocol(prot)
+        output = getattr(prot, outputName, None)
+        self.assertIsNotNone(output, msg)
+        self.assertSetSize(output, 5, msg)
+        self.assertSetSize(output.getImages(), 5236, msg)
+
+    def checkFinal2DClasses2nd(self, prot, outputName="outputClasses", msg=''):
+        prot = self._updateProtocol(prot)
+        output = getattr(prot, outputName, None)
+        self.assertIsNotNone(output, msg)
+        self.assertSetSize(output, 5, msg)
+
+    def checkResults(self, prot, size, msg=''):
+        t0 = time.time()
+        while not prot.isFinished():
+            # Time out 4 minutes, just in case
+            tdelta = time.time() - t0
+            if tdelta > 4 * 60:
+                break
+            prot = self._updateProtocol(prot)
+            time.sleep(2)
+
+        self.assertSetSize(prot.outputClasses, size, msg)
+
 if __name__ == "__main__":
     if len(sys.argv) > 1:
         className = sys.argv[1]
diff --git a/xmipp3/tests/test_protocols_xmipp_3d.py b/xmipp3/tests/test_protocols_xmipp_3d.py
index 55982c2cc..99bf1a7ca 100644
--- a/xmipp3/tests/test_protocols_xmipp_3d.py
+++ b/xmipp3/tests/test_protocols_xmipp_3d.py
@@ -553,6 +553,45 @@ def testSingleFactorAndCrop(self):
         self.assertTrue(outV.equalAttributes(
             inV, ignore=['_index', '_filename', '_samplingRate'], verbose=True))
 
+        outHeader = Ccp4Header(outV.getFileName(), readHeader=True)
+        self.assertAlmostEqual(round(outHeader.getSampling()[0], 3), inV.getSamplingRate() * 2)
+
+    def testSingleResizeSamplingRateAndCrop(self):
+        inV = self.protImport2.outputVolume  # short notation
+        newSR = inV.getSamplingRate() * 2
+        outV = self.launchSingle(doResize=True,
+                                 resizeOption=xrh.RESIZE_SAMPLINGRATE,
+                                 resizeSamplingRate=newSR,
+                                 doWindow=True,
+                                 windowOperation=xrh.WINDOW_OP_CROP)
+
+        self.assertEqual(inV.getDim()[0] * 0.5, outV.getDim()[0])
+        self.assertAlmostEqual(outV.getSamplingRate(), newSR)
+
+        outHeader = Ccp4Header(outV.getFileName(), readHeader=True)
+        self.assertAlmostEqual(round(outHeader.getSampling()[0], 3), newSR)
+        
+        self.assertTrue(outV.equalAttributes(
+            inV, ignore=['_index', '_filename', '_samplingRate'], verbose=True))
+        
+    def testSingleOnlyCropWindow(self):
+        inV = self.protImport2.outputVolume  # short notation
+        newSR = inV.getSamplingRate() # no change
+        newDim = inV.getDim()[0] * 0.5
+        outV = self.launchSingle(doResize=False,
+                                 doWindow=True,
+                                 windowOperation=xrh.WINDOW_OP_WINDOW,
+                                 windowSize=newDim)
+
+        self.assertEqual(newDim, outV.getDim()[0])
+        self.assertAlmostEqual(outV.getSamplingRate(), newSR)
+
+        outHeader = Ccp4Header(outV.getFileName(), readHeader=True)
+        self.assertAlmostEqual(round(outHeader.getSampling()[0], 3), newSR)
+        
+        self.assertTrue(outV.equalAttributes(
+            inV, ignore=['_index', '_filename', '_samplingRate'], verbose=True))
+
     def testSinglePyramid(self):
         inV = self.protImport2.outputVolume  # short notation
         outV = self.launchSingle(doResize=True, resizeOption=xrh.RESIZE_PYRAMID,
@@ -2715,6 +2754,14 @@ def testXmippResolutionBfactor(self):
                                                  fscResolution=8.35)
         self.launchProtocol(protbfactorResolution)
 
+        # Protocol local resolution/local bfactor with assuming already normalized
+        print("Protocol local resolution/local bfactor")
+        protbfactorResolution = self.newProtocol(XmippProtbfactorResolution,
+                                                 normalizeResolution=False,
+                                                 pdbfile=protImportPdb.outputPdb,
+                                                 normalizedMap=protCreateMap.resolution_Volume)
+        self.launchProtocol(protbfactorResolution)
+
 class TestXmippLocalVolAdjust(TestXmippBase):
 
     @classmethod
diff --git a/xmipp3/tests/test_protocols_xmipp_mics.py b/xmipp3/tests/test_protocols_xmipp_mics.py
index 098a65a0c..f7903238f 100644
--- a/xmipp3/tests/test_protocols_xmipp_mics.py
+++ b/xmipp3/tests/test_protocols_xmipp_mics.py
@@ -1185,7 +1185,7 @@ def setUpClass(cls):
                                                  magnification=56000)
 
     def testTiltAnalysis(self):
-        protTilt = XmippProtTiltAnalysis(window_size=512, objective_resolution=8, objLabel=
+        protTilt = XmippProtTiltAnalysis(autoWindow=False, window_size=512, objective_resolution=8, objLabel=
                                         'Tilt analysis window 512 obj resolution 8A')
         protTilt.inputMicrographs.set(self.protImport.outputMicrographs)
         self.proj.launchProtocol(protTilt, wait=True)
@@ -1195,7 +1195,7 @@ def testTiltAnalysis(self):
 
 
     def testTiltAnalysis2(self):
-        protTilt2 = XmippProtTiltAnalysis(window_size=600, objective_resolution=4,
+        protTilt2 = XmippProtTiltAnalysis(autoWindow=False, window_size=600, objective_resolution=4,
                                           meanCorr_threshold=0.9,
                                           objLabel='Tilt analysis window 600 obj resolution 4A and threshold 0.9')
         protTilt2.inputMicrographs.set(self.protImport.outputMicrographs)
diff --git a/xmipp3/tests/test_protocols_xmipp_movies.py b/xmipp3/tests/test_protocols_xmipp_movies.py
index 5b8e81a0e..41146b7a7 100644
--- a/xmipp3/tests/test_protocols_xmipp_movies.py
+++ b/xmipp3/tests/test_protocols_xmipp_movies.py
@@ -572,19 +572,11 @@ def runAlignMovMics(cls):  # do SAVE averaged mics and dose weighted mics
                                     doLocalAlignment=False, useGpu=True,
                                     objLabel='Movie alignment (SAVE mic)',
                                     doSaveAveMic=True,
-                                    binFactor=2)
+                                    binFactor=1)
         protAlign.inputMovies.set(cls.protImport.outputMovies)
         cls.launchProtocol(protAlign)
 
-        protAlign2 = cls.newProtocol(XmippProtOFAlignment,
-                                     doApplyDoseFilter=True,
-                                     alignFrame0=1, alignFrameN=0,
-                                     objLabel='Optical alignment (DW mic)',
-                                     doSaveMovie=True)
-        protAlign2.inputMovies.set(protAlign.outputMovies)
-        cls.launchProtocol(protAlign2)
-
-        return protAlign.outputMovies, protAlign2.outputMovies
+        return protAlign.outputMovies
 
     @classmethod
     def setUpClass(cls):
@@ -593,9 +585,9 @@ def setUpClass(cls):
         fn = 'Falcon_2012_06_12-*0_movie.mrcs'
         cls.protImport = cls.runImportMovies(cls.ds.getFile(fn))
         cls.alignedMovies = cls.runAlignMovies()
-        cls.alignedMovMics, cls.alignedMovDwMics = cls.runAlignMovMics()
+        cls.alignedMovMics = cls.runAlignMovMics()
 
-    def _checkMaxShiftFiltering(self, protocol, label, hasMic, hasDw=False, results=[]):
+    def _checkMaxShiftFiltering(self, protocol, label, hasMic, hasDw=False, noBin=False, results=[]):
         """ Check if outputSets are right.
               If hasMic=True then it's checked that micrographs are generated.
               results = [True, False]  # means first movie should pass whereas second not
@@ -608,7 +600,7 @@ def assertOutput(outputName, ids=[1, 2]):
             targetSamplingRate = protocol.inputMovies.get().getSamplingRate()
             if 'Micrographs' in outputName:
                 # The Mics has n=1
-                if hasDw:
+                if hasDw or noBin:
                     inputDim = (1950, 1950, 1)
                 else:
                     # The corr. align. prot. crops the micrographs and is binned
@@ -702,7 +694,7 @@ def testFilterFrame(self):
         self._checkMaxShiftFiltering(protNoMic, label, hasMic=False, results=[False, True])
 
         protDoMic = self.doFilter(self.alignedMovMics, rejType, label, 0.11) # roughly half the precision of the non-binned version
-        self._checkMaxShiftFiltering(protDoMic, label, hasMic=True, results=[False, True])
+        self._checkMaxShiftFiltering(protDoMic, label, noBin=True, hasMic=True, results=[False, True])
 
     def testFilterMovie(self): 
         """ This must discard the second movie for a Global shift.
@@ -713,8 +705,8 @@ def testFilterMovie(self):
         protNoMic = self.doFilter(self.alignedMovies, rejType, label)
         self._checkMaxShiftFiltering(protNoMic, label, hasMic=False, results=[False, True])
 
-        protDoMic = self.doFilter(self.alignedMovMics, rejType, label, 0.11, 0.16) # roughly half the precision of the non-binned version
-        self._checkMaxShiftFiltering(protDoMic, label, hasMic=True, results=[False, True])
+        protDoMic = self.doFilter(self.alignedMovMics, rejType, label, 0.11, 0.158) # roughly half the precision of the non-binned version
+        self._checkMaxShiftFiltering(protDoMic, label, noBin=True, hasMic=True, results=[True, False])
 
     def testFilterAnd(self): 
         """ This must discard the second movie for AND.
@@ -725,8 +717,8 @@ def testFilterAnd(self):
         protNoMic = self.doFilter(self.alignedMovies, rejType, label) 
         self._checkMaxShiftFiltering(protNoMic, label, hasMic=False, results=[False, True])
 
-        protDoMic = self.doFilter(self.alignedMovMics, rejType, label, 0.11, 0.16) # roughly half the precision of the non-binned version
-        self._checkMaxShiftFiltering(protDoMic, label, hasMic=True, results=[False, True])
+        protDoMic = self.doFilter(self.alignedMovMics, rejType, label, 0.1, 0.158) # roughly half the precision of the non-binned version
+        self._checkMaxShiftFiltering(protDoMic, label, noBin=True, hasMic=True, results=[True, False])
 
     def testFilterOrFrame(self):
         """ This must discard the second movie for OR (Frame).
@@ -738,7 +730,7 @@ def testFilterOrFrame(self):
         self._checkMaxShiftFiltering(protNoMic, label, hasMic=False, results=[False, True])
 
         protDoMic = self.doFilter(self.alignedMovMics, rejType, label, 0.11, 1.0) # roughly half the precision of the non-binned version
-        self._checkMaxShiftFiltering(protDoMic, label, hasMic=True, results=[False, True])
+        self._checkMaxShiftFiltering(protDoMic, label, noBin=True,  hasMic=True, results=[False, True])
 
     def testFilterOrMovie(self): 
         """ This must discard the second movie for OR (Movie).
@@ -749,8 +741,8 @@ def testFilterOrMovie(self):
         protNoMic = self.doFilter(self.alignedMovies, rejType, label, mxFm=1)
         self._checkMaxShiftFiltering(protNoMic, label, hasMic=False, results=[False, True])
 
-        protDoMic = self.doFilter(self.alignedMovMics, rejType, label, 1.0, 0.16) # roughly half the precision of the non-binned version
-        self._checkMaxShiftFiltering(protDoMic, label, hasMic=True, results=[False, True])
+        protDoMic = self.doFilter(self.alignedMovMics, rejType, label, 1.0, 0.158) # roughly half the precision of the non-binned version
+        self._checkMaxShiftFiltering(protDoMic, label, noBin=True, hasMic=True, results=[True, False])
 
     def testFilterOrBoth(self): 
         """ This must discard the second movie for OR (both).
@@ -762,7 +754,7 @@ def testFilterOrBoth(self):
         self._checkMaxShiftFiltering(protNoMic, label, hasMic=False, results=[False, True])
 
         protDoMic = self.doFilter(self.alignedMovMics, rejType, label, 0.11, 0.16) # roughly half the precision of the non-binned version)
-        self._checkMaxShiftFiltering(protDoMic, label, hasMic=True, results=[False, True])
+        self._checkMaxShiftFiltering(protDoMic, label,  noBin=True, hasMic=True, results=[False, True])
 
     def testFilterRejectBoth(self):
         """ This must discard both movies.
@@ -774,7 +766,7 @@ def testFilterRejectBoth(self):
         self._checkMaxShiftFiltering(protNoMic, label, hasMic=False, results=[False, False])
 
         protDoMic = self.doFilter(self.alignedMovMics, rejType, label, mxMo=0.01)
-        self._checkMaxShiftFiltering(protDoMic, label, hasMic=True, results=[False, False])
+        self._checkMaxShiftFiltering(protDoMic, label,  noBin=True, hasMic=True, results=[False, False])
 
     def testFilterAcceptBoth(self):
         """ This must accept both movies.
@@ -786,11 +778,7 @@ def testFilterAcceptBoth(self):
         self._checkMaxShiftFiltering(protNoMic, label, hasMic=False, results=[True, True])
 
         protDoMic = self.doFilter(self.alignedMovMics, rejType, label, mxMo=5)
-        self._checkMaxShiftFiltering(protDoMic, label, hasMic=True, results=[True, True])
-
-        protDoMic = self.doFilter(self.alignedMovDwMics, rejType, label)
-        self._checkMaxShiftFiltering(protDoMic, label, hasMic=True,
-                                     hasDw=True, results=[True, True])
+        self._checkMaxShiftFiltering(protDoMic, label,  noBin=True, hasMic=True, results=[True, True])
 
 
 class TestMovieDoseAnalysis(BaseTest):
@@ -869,4 +857,113 @@ def testDoseAnalysisFiltering(self):
 
         sizeDiscarded = protPoisson.outputMoviesDiscarded.getSize()
         self.assertEqual(sizeDiscarded, 5, 'Number of accepted movies must be 5 and its '
-                                            '%d' % sizeDiscarded)
\ No newline at end of file
+                                            '%d' % sizeDiscarded)
+
+
+class TestMovieAlignmentConsensus(BaseTest):
+    @classmethod
+    def setData(cls):
+        cls.ds = DataSet.getDataSet('movies')
+
+    @classmethod
+    def runImportMovies(cls, pattern, **kwargs):
+        """ Run an Import movies protocol. """
+        # We have two options: passes the SamplingRate or
+        # the ScannedPixelSize + microscope magnification
+        params = {'samplingRate': 1.14,
+                  'voltage': 300,
+                  'sphericalAberration': 2.7,
+                  'magnification': 50000,
+                  'scannedPixelSize': None,
+                  'filesPattern': pattern,
+                  'dosePerFrame': 123
+                  }
+        if 'samplingRate' not in kwargs:
+            del params['samplingRate']
+            params['samplingRateMode'] = 0
+        else:
+            params['samplingRateMode'] = 1
+
+        params.update(kwargs)
+
+        protImport = cls.newProtocol(ProtImportMovies, **params)
+        cls.launchProtocol(protImport)
+        return protImport
+
+    @classmethod
+    def runAlignMovies1(cls):
+        protAlign = cls.newProtocol(XmippProtFlexAlign,
+                                    alignFrame0=1, alignFrameN=0,
+                                    doLocalAlignment=False, useGpu=False,
+                                    objLabel='Reference movie alignment',
+                                    doSaveAveMic=True)
+
+        protAlign.inputMovies.set(cls.protImport.outputMovies)
+        cls.launchProtocol(protAlign)
+        return protAlign
+
+    @classmethod
+    def runAlignMovies2(cls):
+        protAlign2 = cls.newProtocol(XmippProtFlexAlign,
+                                    alignFrame0=1, alignFrameN=0,
+                                    maxResForCorrelation=20,
+                                    doLocalAlignment=False, useGpu=False,
+                                    objLabel='Target movie alignment',
+                                    doSaveAveMic=True)
+
+        protAlign2.inputMovies.set(cls.protImport.outputMovies)
+        cls.launchProtocol(protAlign2)
+        return protAlign2
+
+    @classmethod
+    def setUpClass(cls):
+        setupTestProject(cls)
+        cls.setData()
+        fn = 'Falcon_2012_06_12-*0_movie.mrcs'
+        cls.protImport = cls.runImportMovies(cls.ds.getFile(fn))
+        cls.align1 = cls.runAlignMovies1()
+        cls.align2 = cls.runAlignMovies2()
+
+    def testMovieAlignmentConsensusFiltering1(self):
+        """ This must discard movies by movie alignment consensus.
+        """
+        label = 'Alignment consensus 0.5 correlation limit'
+        protConsensus1 = self.newProtocol(XmippProtConsensusMovieAlignment,
+                                          objLabel=label,
+                                          minConsCorrelation=0.5,
+                                          trajectoryPlot=True
+                                          )
+
+        protConsensus1.inputMovies1.set(self.align1)
+        protConsensus1.inputMovies1.setExtended("outputMovies")
+
+        protConsensus1.inputMovies2.set(self.align2)
+        protConsensus1.inputMovies2.setExtended("outputMovies")
+
+        self.launchProtocol(protConsensus1)
+
+        sizeAccepted = protConsensus1.outputMovies.getSize()
+        self.assertEqual(sizeAccepted, 2, 'Number of accepted movies must be 2 and its '
+                                           '%d' % sizeAccepted)
+
+    def testMovieAlignmentConsensusFiltering2(self):
+        """ This must discard movies by movie alignment consensus.
+        """
+        label = 'Alignment consensus 0.9 correlation limit'
+        protConsensus2 = self.newProtocol(XmippProtConsensusMovieAlignment,
+                                          objLabel=label,
+                                          minConsCorrelation=0.9,
+                                          trajectoryPlot=True
+                                          )
+
+        protConsensus2.inputMovies1.set(self.align1)
+        protConsensus2.inputMovies1.setExtended("outputMovies")
+
+        protConsensus2.inputMovies2.set(self.align2)
+        protConsensus2.inputMovies2.setExtended("outputMovies")
+
+        self.launchProtocol(protConsensus2)
+
+        sizeDiscarded = protConsensus2.outputMoviesDiscarded.getSize()
+        self.assertEqual(sizeDiscarded, 2, 'Number of accepted movies must be 2 and its '
+                                           '%d' % sizeDiscarded)
\ No newline at end of file