SNSPowderReduction.py 75.5 KB
Newer Older
1
2
3
# Mantid Repository : https://github.com/mantidproject/mantid
#
# Copyright © 2018 ISIS Rutherford Appleton Laboratory UKRI,
4
5
#   NScD Oak Ridge National Laboratory, European Spallation Source,
#   Institut Laue - Langevin & CSNS, Institute of High Energy Physics, CAS
6
# SPDX - License - Identifier: GPL - 3.0 +
7
#pylint: disable=invalid-name,no-init,too-many-lines
8
import os
9
from pathlib import Path
10
import mantid.simpleapi as api
11
from mantid.api import mtd, AlgorithmFactory, AnalysisDataService, DistributedDataProcessorAlgorithm, \
Peterson, Peter's avatar
Peterson, Peter committed
12
    FileAction, FileProperty, ITableWorkspaceProperty, MultipleFileProperty, PropertyMode, WorkspaceProperty, \
13
    ITableWorkspace, MatrixWorkspace
14
15
16
from mantid.kernel import (
    ConfigService, Direction, EnabledWhenProperty, FloatArrayProperty, FloatBoundedValidator, IntArrayBoundedValidator,
    IntArrayProperty, Property, PropertyCriterion, PropertyManagerDataService, StringListValidator)
Zhou, Wenduo's avatar
Zhou, Wenduo committed
17
from mantid.dataobjects import SplittersWorkspace  # SplittersWorkspace
18
from mantid.utils import absorptioncorrutils
19
if AlgorithmFactory.exists('GatherWorkspaces'):
20
    HAVE_MPI = True
21
    from mpi4py import MPI
22
    mpiRank = MPI.COMM_WORLD.Get_rank()
23
24
else:
    HAVE_MPI = False
25
    mpiRank = 0 # simplify if clauses
26

27

28
EVENT_WORKSPACE_ID = "EventWorkspace"
29
EXTENSIONS_NXS = ["_event.nxs", ".nxs.h5"]
30

Whitfield, Ross's avatar
Whitfield, Ross committed
31

32
def noRunSpecified(runs):
33
34
35
36
37
38
39
    """ Check whether any run is specified
    Purpose:
        Check whether any run, which is larger than 0, is specified in the input numpy array
    Requirements:
        Input is a numpy array
    Guarantees:
        A boolean is returned
40
41
42
    :param runs:
    :return:
    """
43
    # return True if runs is of size zero
44
    if len(runs) <= 0:
45
        return True
46
47

    # return True if there is one and only one run in runs and it is not greater than 0.
48
49
50
51
52
    if len(runs) == 1:
        # break off the instrument part
        value = int(runs[0].split('_')[-1])
        # -1 turns off the runnumber
        return value <= 0
53

54
55
    return False

Whitfield, Ross's avatar
Whitfield, Ross committed
56

57
58
59
60
61
62
63
64
65
66
67
def get_workspace(workspace_name):
    """
    Purpose: Get the reference of a workspace
    Requirements:
    1. workspace_name is a string
    Guarantees:
    The reference of the workspace with same is returned.
    :exception RuntimeError: a RuntimeError from Mantid will be thrown
    :param workspace_name:
    :return:
    """
68
69
    assert isinstance(workspace_name, str), 'Input workspace name {0} must be a string but not a {1}.' \
                                            ''.format(workspace_name, type(workspace_name))
70
71
72

    return AnalysisDataService.retrieve(workspace_name)

Whitfield, Ross's avatar
Whitfield, Ross committed
73

74
75
76
77
78
79
def is_event_workspace(workspace):
    if isinstance(workspace, str):
        return get_workspace(workspace).id() == EVENT_WORKSPACE_ID
    else:
        return workspace.id() == EVENT_WORKSPACE_ID

Whitfield, Ross's avatar
Whitfield, Ross committed
80

81
def allEventWorkspaces(*args):
82
83
84
85
86
87
88
89
90
91
    """
    Purpose:
        Check whether all the inputs are event workspace
    Requirements:
        Each element in the args is a string as name of workspace
    Guarantees:
        return boolean
    @param args:
    @return:
    """
92
93
    result = True

94
    args = [is_event_workspace(arg) for arg in args]
95
    for arg in args:
96
        result = result and arg
97
98
99

    return result

Whitfield, Ross's avatar
Whitfield, Ross committed
100

101
def getBasename(filename):
102
103
    if type(filename) == list:
        filename = filename[0]
104
105
106
107
108
    name = os.path.split(filename)[-1]
    for extension in EXTENSIONS_NXS:
        name = name.replace(extension, '')
    return name

Savici, Andrei T's avatar
Savici, Andrei T committed
109
#pylint: disable=too-many-instance-attributes
Whitfield, Ross's avatar
Whitfield, Ross committed
110
111


112
class SNSPowderReduction(DistributedDataProcessorAlgorithm):
113
    COMPRESS_TOL_TOF = .01
114
115
116
117
118
119
120
121
122
    _resampleX = None
    _binning = None
    _bin_in_dspace = None
    _instrument = None
    _filterBadPulses = None
    _removePromptPulseWidth = None
    _LRef = None
    _DIFCref = None
    _wavelengthMin = None
123
    _wavelengthMax = None
124
125
    _vanPeakFWHM = None
    _vanSmoothing = None
126
    _vanRadius = None
127
128
129
130
131
    _scaleFactor = None
    _outDir = None
    _outPrefix = None
    _outTypes = None
    _chunks = None
Zhou, Wenduo's avatar
Zhou, Wenduo committed
132
    _splittersWS = None
133
134
135
136
137
138
139
    _splitinfotablews = None
    _normalisebycurrent = None
    _lowResTOFoffset = None
    _focusPos = None
    _charTable = None
    iparmFile = None
    _info = None
140
    _absMethod = None
141
142
    _sampleFormula = None
    _massDensity = None
143
    _containerShape = None
144

145
    def category(self):
146
        return "Diffraction\\Reduction"
147

Nick Draper's avatar
Nick Draper committed
148
149
150
    def seeAlso(self):
        return [ "DiffractionFocussing","AlignAndFocusPowder" ]

151
152
153
    def name(self):
        return "SNSPowderReduction"

154
    def summary(self):
155
156
157
158
        """
        summary of the algorithm
        :return:
        """
Alex Buts's avatar
Alex Buts committed
159
        return "The algorithm used for reduction of powder diffraction data obtained on SNS instruments (e.g. PG3) "
160

161
    def PyInit(self):
162
163
        self.copyProperties('AlignAndFocusPowderFromFiles', ['Filename', 'PreserveEvents'])

164
        self.declareProperty("Sum", False,
165
                             "Sum the runs. Does nothing for characterization runs")
166
        self.declareProperty("PushDataPositive", "None",
167
168
                             StringListValidator(["None", "ResetToZero", "AddMinimum"]),
                             "Add a constant to the data that makes it positive over the whole range.")
169
170
        arrvalidator = IntArrayBoundedValidator(lower=-1)
        self.declareProperty(IntArrayProperty("BackgroundNumber", values=[0], validator=arrvalidator),
171
                             doc="If specified overrides value in CharacterizationRunsFile If -1 turns off correction.")
172
173
174
175
        self.declareProperty(IntArrayProperty("VanadiumNumber", values=[0], validator=arrvalidator),
                             doc="If specified overrides value in CharacterizationRunsFile. If -1 turns off correction.")
        self.declareProperty(IntArrayProperty("VanadiumBackgroundNumber", values=[0], validator=arrvalidator),
                             doc="If specified overrides value in CharacterizationRunsFile. If -1 turns off correction.")
176
        self.declareProperty(FileProperty(name="CalibrationFile",defaultValue="",action=FileAction.OptionalLoad,
177
                                          extensions=[".h5", ".hd5", ".hdf", ".cal"]))  # CalFileName
Whitfield, Ross's avatar
Whitfield, Ross committed
178
        self.declareProperty(FileProperty(name="GroupingFile",defaultValue="",action=FileAction.OptionalLoad,
179
                                          extensions=[".xml"]), "Overrides grouping from CalibrationFile")
Peterson, Peter's avatar
Peterson, Peter committed
180
        self.declareProperty(MultipleFileProperty(name="CharacterizationRunsFile",
Peterson, Peter's avatar
Peterson, Peter committed
181
182
                                                  action=FileAction.OptionalLoad,
                                                  extensions=["txt"]), "File with characterization runs denoted")
183
184
        self.declareProperty(FileProperty(name="ExpIniFilename", defaultValue="", action=FileAction.OptionalLoad,
                                          extensions=[".ini"]))
185
        self.copyProperties('AlignAndFocusPowderFromFiles',
186
187
                            ['LorentzCorrection', 'UnwrapRef', 'LowResRef', 'CropWavelengthMin', 'CropWavelengthMax',
                             'RemovePromptPulseWidth', 'MaxChunkSize'])
188
189
        self.declareProperty(FloatArrayProperty("Binning", values=[0., 0., 0.],
                                                direction=Direction.Input),
190
191
                             "Positive is linear bins, negative is logorithmic")  # Params
        self.copyProperties('AlignAndFocusPowderFromFiles', ['ResampleX'])
Peterson, Peter's avatar
Peterson, Peter committed
192
        self.declareProperty("BinInDspace", True,
193
                             "If all three bin parameters a specified, whether they are in dspace (true) or "
194
                             "time-of-flight (false)")  # DSpacing
195
        # section of vanadium run processing
196
        self.declareProperty("StripVanadiumPeaks", True,
197
198
                             "Subtract fitted vanadium peaks from the known positions.")
        self.declareProperty("VanadiumFWHM", 7, "Default=7")
199
        self.declareProperty("VanadiumPeakTol", 0.05,
Whitfield, Ross's avatar
Whitfield, Ross committed
200
                             "How far from the ideal position a vanadium peak can be during StripVanadiumPeaks. "
Savici, Andrei T's avatar
Savici, Andrei T committed
201
                             "Default=0.05, negative turns off")
202
        self.declareProperty("VanadiumSmoothParams", "20,2", "Default=20,2")
203
        self.declareProperty("VanadiumRadius", .3175, "Radius for CarpenterSampleCorrection")
204
        self.declareProperty("BackgroundSmoothParams", "", "Default=off, suggested 20,2")
205
206

        # filtering
207
        self.declareProperty("FilterBadPulses", 95.,  # different default value
208
                             doc="Filter out events measured while proton charge is more than 5% below average")
209
        self.declareProperty("ScaleData", defaultValue=1., validator=FloatBoundedValidator(lower=0., exclusive=True),
210
211
                             doc="Constant to multiply the data before writing out. This does not apply to "
                                 "PDFgetN files.")
212
213
214
        self.declareProperty("OffsetData", defaultValue=0., validator=FloatBoundedValidator(lower=0., exclusive=False),
                             doc="Constant to add to the data before writing out. This does not apply to "
                                 "PDFgetN files.")
215
        self.declareProperty("SaveAs", "gsas",
216
217
                             "List of all output file types. Allowed values are 'fullprof', 'gsas', 'nexus', "
                             "'pdfgetn', and 'topas'")
218
        self.declareProperty("OutputFilePrefix", "", "Overrides the default filename for the output file (Optional).")
219
220
        self.declareProperty(FileProperty(name="OutputDirectory", defaultValue="",action=FileAction.Directory))

221
        # Caching options
222
        self.copyProperties('AlignAndFocusPowderFromFiles', 'CacheDir')
223
224
225
        self.declareProperty('CleanCache', False, 'Remove all cache files within CacheDir')
        self.setPropertySettings('CleanCache', EnabledWhenProperty('CacheDir', PropertyCriterion.IsNotDefault))
        property_names = ('CacheDir', 'CleanCache')
226
        [self.setPropertyGroup(name, 'Caching') for name in property_names]
227

228
        self.declareProperty("FinalDataUnits", "dSpacing", StringListValidator(["dSpacing","MomentumTransfer"]))
229

230
231
        # absorption correction
        self.declareProperty("TypeOfCorrection", "None",
232
233
234
                             StringListValidator(["None", "SampleOnly", "SampleAndContainer", "FullPaalmanPings"]),
                             doc="Specifies the Absorption Correction terms to calculate, if any.")
        self.declareProperty("SampleFormula", "", doc="Chemical formula of the sample")
Kendrick, Coleman's avatar
Kendrick, Coleman committed
235
        self.declareProperty("MeasuredMassDensity", defaultValue=0.1,
236
237
                             validator=FloatBoundedValidator(lower=0., exclusive=True),
                             doc="Measured mass density of sample in g/cc")  # in g/cc, way to validate?
238
239
        self.declareProperty("SampleNumberDensity", defaultValue=Property.EMPTY_DBL,
                             doc="Number density of the sample in number of atoms per cubic Angstrom will be used instead of calculated")
240
        self.declareProperty("ContainerShape", defaultValue="PAC06", doc="Defines the container geometry")
241
242
243
244
245
246
        self.declareProperty("ContainerScaleFactor", defaultValue=1.0,
                             validator=FloatBoundedValidator(lower=0),
                             doc="Factor to scale the container data")
        self.copyProperties("AbsorptionCorrection", "ElementSize")
        self.declareProperty("NumWavelengthBins", defaultValue=1000,
                             doc="Number of wavelength bin to calculate the for absorption correction")
247

248
249
250
251
252
253
        workspace_prop = WorkspaceProperty('SplittersWorkspace', '', Direction.Input, PropertyMode.Optional)
        self.declareProperty(workspace_prop, "Splitters workspace for split event workspace.")
        # replaced for matrix workspace, SplittersWorkspace and table workspace
        # tableprop = ITableWorkspaceProperty("SplittersWorkspace", "", Direction.Input, PropertyMode.Optional)
        # self.declareProperty(tableprop, "Splitters workspace for split event workspace.")

Zhou, Wenduo's avatar
Zhou, Wenduo committed
254
255
        infotableprop = ITableWorkspaceProperty("SplitInformationWorkspace", "", Direction.Input, PropertyMode.Optional)
        self.declareProperty(infotableprop, "Name of table workspace containing information for splitters.")
256

257
        self.declareProperty("LowResolutionSpectraOffset", -1,
258
                             "If larger and equal to 0, then process low resolution TOF and offset is the spectra "
259
                             "number. Otherwise, ignored.")  # LowResolutionSpectraOffset
260

261
262
        self.declareProperty("NormalizeByCurrent", True, "Normalize by current")

263
264
        self.declareProperty("CompressTOFTolerance", 0.01, "Tolerance to compress events in TOF.")

265
        self.copyProperties('AlignAndFocusPowderFromFiles', ['FrequencyLogNames', 'WaveLengthLogNames'])
266

267
        return
268

269
270
271
    def validateInputs(self):
        issues = dict()

Kendrick, Coleman's avatar
Kendrick, Coleman committed
272
273
274
275
276
        # If doing absorption correction, make sure the sample formula is correct
        if self.getProperty("TypeOfCorrection").value != "None":
            if self.getProperty("SampleFormula").value == '':
                issues['SampleFormula'] = "A sample formula must be provided."

277
        # The provided cache directory does not exist
278
        cache_dir = self.getProperty('CacheDir').value  # absolute or relative path, as a string
Jose Borreguero's avatar
Jose Borreguero committed
279
        if bool(cache_dir) and Path(cache_dir).exists() is False:
280
            issues['CacheDir'] = f'Directory {cache_dir} does not exist'
281
282

        # We cannot clear the cache if property "CacheDir" has not been set
283
284
        if self.getProperty('CleanCache').value and not bool(self.getProperty('CacheDir').value):
            issues['CleanCache'] = f'Property "CacheDir" must be set in order to clean the cache'
285

286
287
        return issues

288
    #pylint: disable=too-many-locals,too-many-branches,too-many-statements
289
    def PyExec(self):  # noqa
290
291
        """ Main execution body
        """
292
        # get generic information
293
        self._loadCharacterizations()
294
295
296
297
298
299
300
301
302
303
        self._resampleX = self.getProperty("ResampleX").value
        if self._resampleX != 0.:
            self._binning = [0.]
        else:
            self._binning = self.getProperty("Binning").value
            if len(self._binning) != 1 and len(self._binning) != 3:
                raise RuntimeError("Can only specify (width) or (start,width,stop) for binning. Found %d values." % len(self._binning))
            if len(self._binning) == 3:
                if self._binning[0] == 0. and self._binning[1] == 0. and self._binning[2] == 0.:
                    raise RuntimeError("Failed to specify the binning")
304
305
306
        self._bin_in_dspace = self.getProperty("BinInDspace").value
        self._filterBadPulses = self.getProperty("FilterBadPulses").value
        self._removePromptPulseWidth = self.getProperty("RemovePromptPulseWidth").value
307
        self._lorentz = self.getProperty("LorentzCorrection").value
308
309
310
        self._LRef = self.getProperty("UnwrapRef").value
        self._DIFCref = self.getProperty("LowResRef").value
        self._wavelengthMin = self.getProperty("CropWavelengthMin").value
311
        self._wavelengthMax = self.getProperty("CropWavelengthMax").value
312
313
        self._vanPeakFWHM = self.getProperty("VanadiumFWHM").value
        self._vanSmoothing = self.getProperty("VanadiumSmoothParams").value
314
        self._vanRadius = self.getProperty("VanadiumRadius").value
315
        self.calib = self.getProperty("CalibrationFile").value
316
        self._scaleFactor = self.getProperty("ScaleData").value
317
        self._offsetFactor = self.getProperty("OffsetData").value
318
        self._outDir = self.getProperty("OutputDirectory").value
319
        # Caching options
320
321
322
        self._cache_dir = self.getProperty("CacheDir").value
        self._clean_cache = self.getProperty("CleanCache").value

Peterson, Peter's avatar
Peterson, Peter committed
323
        self._outPrefix = self.getProperty("OutputFilePrefix").value.strip()
324
        self._outTypes = self.getProperty("SaveAs").value.lower()
325
326
327
        self._absMethod = self.getProperty("TypeOfCorrection").value
        self._sampleFormula = self.getProperty("SampleFormula").value
        self._massDensity = self.getProperty("MeasuredMassDensity").value
328
        self._numberDensity = self.getProperty("SampleNumberDensity").value
Whitfield, Ross's avatar
Whitfield, Ross committed
329
        self._containerShape = self.getProperty("ContainerShape").value
330
331
332
        self._containerScaleFactor = self.getProperty("ContainerScaleFactor").value
        self._elementSize = self.getProperty("ElementSize").value
        self._num_wl_bins = self.getProperty("NumWavelengthBins").value
333

334
        samRuns = self._getLinearizedFilenames("Filename")
335
336
        self._determineInstrument(samRuns[0])

337
        preserveEvents = self.getProperty("PreserveEvents").value
338
        if HAVE_MPI and preserveEvents:
339
            self.log().warning("preserveEvents set to False for MPI tasks.")
340
            preserveEvents = False
341
        self._info = None
342
        self._chunks = self.getProperty("MaxChunkSize").value
343

Zhou, Wenduo's avatar
Zhou, Wenduo committed
344
345
346
347
348
        # define splitters workspace and filter wall time
        self._splittersWS = self.getProperty("SplittersWorkspace").value
        if self._splittersWS is not None:
            # user specifies splitters workspace
            self.log().information("SplittersWorkspace is %s" % (str(self._splittersWS)))
349
            if len(samRuns) != 1:
Zhou, Wenduo's avatar
Zhou, Wenduo committed
350
351
352
                # TODO/FUTURE - This should be supported
                raise RuntimeError("Reducing data with splitters cannot happen when there are more than 1 sample run.")
            # define range of wall-time to import data
353
            sample_time_filter_wall = self._get_time_filter_wall(self._splittersWS.name(), samRuns[0])
354
            self.log().information("The time filter wall is %s" %(str(sample_time_filter_wall)))
355
        else:
356
            sample_time_filter_wall = (0.0, 0.0)
357
            self.log().information("SplittersWorkspace is None, and thus there is NO time filter wall. ")
358

359
        self._splitinfotablews = self.getProperty("SplitInformationWorkspace").value
Zhou, Wenduo's avatar
Zhou, Wenduo committed
360

361
362
        self._normalisebycurrent = self.getProperty("NormalizeByCurrent").value

363
        # Tolerance for compress TOF event.  If given a negative value, then use default 0.01
364
        self.COMPRESS_TOL_TOF = float(self.getProperty("CompressTOFTolerance").value)
365
        if self.COMPRESS_TOL_TOF < -0.:
366
            self.COMPRESS_TOL_TOF = 0.01
367

368
369
370
371
        # Clean the cache directory if so requested
        if self._clean_cache:
            api.CleanFileCache(CacheDir=self._cache_dir, AgeInDays=0)

372
        # Process data
Zhou, Wenduo's avatar
Zhou, Wenduo committed
373
374
        # List stores the workspacename of all data workspaces that will be converted to d-spacing in the end.
        workspacelist = []
375
        samwksplist = []
376

377
        self._lowResTOFoffset = self.getProperty("LowResolutionSpectraOffset").value
378
        focuspos = self._focusPos
379
380
        if self._lowResTOFoffset >= 0:
            # Dealing with the parameters for editing instrument parameters
Peterson, Peter's avatar
Peterson, Peter committed
381
            if "PrimaryFlightPath" in focuspos:
382
383
384
385
386
387
388
                l1 = focuspos["PrimaryFlightPath"]
                if l1 > 0:
                    specids = focuspos['SpectrumIDs'][:]
                    l2s = focuspos['L2'][:]
                    polars = focuspos['Polar'][:]
                    phis = focuspos['Azimuthal'][:]

389
390
                    specids.extend(specids)
                    l2s.extend(l2s)
391
392
393
394
395
396
397
398
399
                    polars.extend(polars)
                    phis.extend(phis)

                    focuspos['SpectrumIDs'] = specids
                    focuspos['L2'] = l2s
                    focuspos['Polar'] = polars
                    focuspos['Azimuthal'] = phis
        # ENDIF

400
        # calculate absorption from first sample run
401
402
        metaws = None
        if self._absMethod != "None" and self._info is None:
403
404
405
            absName = '__{}_abs'.format(getBasename(samRuns[0]))
            api.Load(Filename=samRuns[0], OutputWorkspace=absName, MetaDataOnly=True)
            self._info = self._getinfo(absName)
406
            metaws = absName
407
408
409
410
411
412
413
414
415
416
417
418
419
        # NOTE: inconsistent naming among different methods
        #       -> adding more comments to help clarify
        a_sample, a_container = absorptioncorrutils.calculate_absorption_correction(
            samRuns[0],  # filename: File to be used for absorption correction
            self._absMethod,  # [None, SampleOnly, SampleAndContainer, FullPaalmanPings]
            self._info,  # PropertyManager of run characterizations
            self._sampleFormula,  # Material for absorption correction
            self._massDensity,  # Mass density of the sample
            self._numberDensity,  # Optional number density of sample to be added
            self._containerShape,  # Shape definition of container
            self._num_wl_bins,  # Number of bins: len(ws.readX(0))-1
            self._elementSize,  # Size of one side of the integration element cube in mm
            metaws,  # Optional workspace containing metadata
Zhang, Chen's avatar
Zhang, Chen committed
420
            self.getProperty("CacheDir").value,  # Cache dir for absoption correction workspace
421
            "SHA", # Use cache files named with sha rather than a filename prefix
422
        )
423

424
425
        if self.getProperty("Sum").value and len(samRuns) > 1:
            self.log().information('Ignoring value of "Sum" property')
426
            # Sum input sample runs and then do reduction
Zhou, Wenduo's avatar
Zhou, Wenduo committed
427
            if self._splittersWS is not None:
428
429
                raise NotImplementedError("Summing spectra and filtering events are not supported simultaneously.")

430
            sam_ws_name = self._focusAndSum(samRuns, preserveEvents=preserveEvents, absorptionWksp=a_sample)
Zhou, Wenduo's avatar
Zhou, Wenduo committed
431
432
            assert isinstance(sam_ws_name, str), 'Returned from _focusAndSum() must be a string but not' \
                                                 '%s. ' % str(type(sam_ws_name))
433

Zhou, Wenduo's avatar
Zhou, Wenduo committed
434
435
            workspacelist.append(sam_ws_name)
            samwksplist.append(sam_ws_name)
436
        # ENDIF (SUM)
437
438
439
440
        else:
            # Process each sample run
            for sam_run_number in samRuns:
                # first round of processing the sample
441
                self._info = None
442
443
                if sample_time_filter_wall[0] == 0. and sample_time_filter_wall[-1] == 0. \
                        and self._splittersWS is None:
444
                    returned = self._focusAndSum([sam_run_number], preserveEvents=preserveEvents, absorptionWksp=a_sample)
445
446
447
448
                else:
                    returned = self._focusChunks(sam_run_number, sample_time_filter_wall,
                                                 splitwksp=self._splittersWS,
                                                 preserveEvents=preserveEvents)
449

450
                if isinstance(returned, list):
451
452
                    # Returned with a list of workspaces
                    focusedwksplist = returned
Zhou, Wenduo's avatar
Zhou, Wenduo committed
453
454
455
456
457
                    for sam_ws_name in focusedwksplist:
                        assert isinstance(sam_ws_name, str), 'Impossible to have a non-string value in ' \
                                                             'returned focused workspaces\' names.'
                        samwksplist.append(sam_ws_name)
                        workspacelist.append(sam_ws_name)
458
                    # END-FOR
459
                else:
Zhou, Wenduo's avatar
Zhou, Wenduo committed
460
461
462
463
464
                    # returned as a single workspace
                    sam_ws_name = returned
                    assert isinstance(sam_ws_name, str)
                    samwksplist.append(sam_ws_name)
                    workspacelist.append(sam_ws_name)
465
                # ENDIF
466
467
            # END-FOR
        # ENDIF (Sum data or not)
468

469
        for (samRunIndex, sam_ws_name) in enumerate(samwksplist):
Whitfield, Ross's avatar
Whitfield, Ross committed
470
            assert isinstance(sam_ws_name, str), 'Assuming that samRun is a string. But it is %s' % str(type(sam_ws_name))
471
            if is_event_workspace(sam_ws_name):
472
                self.log().information('Sample Run %s:  starting number of events = %d.' % (
473
                    sam_ws_name, get_workspace(sam_ws_name).getNumberEvents()))
474

475
476
            # Get run information
            self._info = self._getinfo(sam_ws_name)
477
478

            # process the container
479
            can_run_numbers = self._info["container"].value
480
            can_run_numbers = ['%s_%d' % (self._instrument, value) for value in can_run_numbers]
481
            can_run_ws_name = self._process_container_runs(can_run_numbers, samRunIndex, preserveEvents, absorptionWksp=a_container)
482
            if can_run_ws_name is not None:
483
                workspacelist.append(can_run_ws_name)
484
485

            # process the vanadium run
486
            van_run_number_list = self._info["vanadium"].value
487
            van_run_number_list = ['%s_%d' % (self._instrument, value) for value in van_run_number_list]
488
489
            van_specified = not noRunSpecified(van_run_number_list)
            if van_specified:
490
                van_run_ws_name = self._process_vanadium_runs(van_run_number_list, samRunIndex)
491
                workspacelist.append(van_run_ws_name)
492
            else:
493
                van_run_ws_name = None
494

495
            # return if MPI is used and there is more than 1 processor
496
497
            if mpiRank > 0:
                return
498
499

            # return if there is no sample run
500
501
502
503
            # Note: sample run must exist in logic
            # VZ: Discuss with Pete
            # if sam_ws_name == 0:
            #   return
504
505

            # the final bit of math to remove container run and vanadium run
506
            if can_run_ws_name is not None and self._containerScaleFactor != 0:
507
                # must convert the sample to a matrix workspace if the can run isn't one
508
                if not allEventWorkspaces(can_run_ws_name, sam_ws_name):
509
510
                    api.ConvertToMatrixWorkspace(InputWorkspace=sam_ws_name,
                                                 OutputWorkspace=sam_ws_name)
511

512
                # remove container run
513
                api.RebinToWorkspace(WorkspaceToRebin=can_run_ws_name,
Peterson, Peter's avatar
Peterson, Peter committed
514
515
                                     WorkspaceToMatch=sam_ws_name,
                                     OutputWorkspace=can_run_ws_name)
516
517
518
                api.Scale(InputWorkspace=can_run_ws_name,
                          OutputWorkspace=can_run_ws_name,
                          Factor=self._containerScaleFactor)
519
520
521
                api.Minus(LHSWorkspace=sam_ws_name,
                          RHSWorkspace=can_run_ws_name,
                          OutputWorkspace=sam_ws_name)
522
                # compress event if the sample run workspace is EventWorkspace
523
                if is_event_workspace(sam_ws_name) and self.COMPRESS_TOL_TOF > 0.:
524
525
526
                    api.CompressEvents(InputWorkspace=sam_ws_name,
                                       OutputWorkspace=sam_ws_name,
                                       Tolerance=self.COMPRESS_TOL_TOF)  # 10ns
527
528
529
530
                # canRun = str(canRun)

            if van_run_ws_name is not None:
                # subtract vanadium run from sample run by division
531
532
                num_hist_sam = get_workspace(sam_ws_name).getNumberHistograms()
                num_hist_van = get_workspace(van_run_ws_name).getNumberHistograms()
533
534
                assert num_hist_sam == num_hist_van, \
                    'Number of histograms of sample %d is not equal to van %d.' % (num_hist_sam, num_hist_van)
535
536
537
                api.Divide(LHSWorkspace=sam_ws_name,
                           RHSWorkspace=van_run_ws_name,
                           OutputWorkspace=sam_ws_name)
538
                normalized = True
539
540
                van_run_ws = get_workspace(van_run_ws_name)
                sam_ws = get_workspace(sam_ws_name)
541
                sam_ws.getRun()['van_number'] = van_run_ws.getRun()['run_number'].value
542
                # van_run_number = str(van_run_number)
543
544
            else:
                normalized = False
545

546
            # Compress the event again
547
            if is_event_workspace(sam_ws_name) and self.COMPRESS_TOL_TOF > 0.:
548
549
550
                api.CompressEvents(InputWorkspace=sam_ws_name,
                                   OutputWorkspace=sam_ws_name,
                                   Tolerance=self.COMPRESS_TOL_TOF)  # 5ns/
551
552

            # write out the files
553
            # FIXME - need documentation for mpiRank
554
            if mpiRank == 0:
555
                if self._scaleFactor != 1.:
556
                    api.Scale(sam_ws_name, Factor=self._scaleFactor, OutputWorkspace=sam_ws_name)
557
558
559
560
561
562
563
564
565
566
567
568
569
                if self._offsetFactor != 0.:
                    api.ConvertToMatrixWorkspace(InputWorkspace=sam_ws_name,
                                                 OutputWorkspace=sam_ws_name)
                    api.Scale(sam_ws_name, Factor=self._offsetFactor, OutputWorkspace=sam_ws_name,
                              Operation='Add')
                # make sure there are no negative values - gsas hates them
                if self.getProperty("PushDataPositive").value != "None":
                    addMin = self.getProperty("PushDataPositive").value == "AddMinimum"
                    api.ResetNegatives(InputWorkspace=sam_ws_name,
                                       OutputWorkspace=sam_ws_name,
                                       AddMinimum=addMin,
                                       ResetValue=0.)

570
                self._save(sam_ws_name, self._info, normalized, False)
571
        # ENDFOR
572

573
        # convert everything into d-spacing as the final units
574
        if mpiRank == 0:
575
            workspacelist = set(workspacelist) # only do each workspace once
576
            for wksp in workspacelist:
577
578
579
                api.ConvertUnits(InputWorkspace=wksp, OutputWorkspace=wksp,
                                 Target=self.getProperty("FinalDataUnits").value)

580
581
582
583
                propertyName = "OutputWorkspace%s" % wksp
                if not self.existsProperty(propertyName):
                    self.declareProperty(WorkspaceProperty(propertyName, wksp, Direction.Output))
                self.setProperty(propertyName, wksp)
584

585
586
        return

587
    def _getLinearizedFilenames(self, propertyName):
Pete Peterson's avatar
Pete Peterson committed
588
        runnumbers = self.getProperty(propertyName).value
589
590
591
592
593
594
595
596
        linearizedRuns = []
        for item in runnumbers:
            if type(item) == list:
                linearizedRuns.extend(item)
            else:
                linearizedRuns.append(item)
        return linearizedRuns

597
598
    def _determineInstrument(self, filename):
        name = getBasename(filename)
599
600
        parts = name.split('_')
        self._instrument = ConfigService.getInstrument(parts[0]).shortName()  # only works for instruments without '_'
601

602
    def _loadCharacterizations(self):
603
        self._focusPos = {}
604
        charFilename = self.getProperty("CharacterizationRunsFile").value
Peterson, Peter's avatar
Peterson, Peter committed
605
        charFilename = ','.join(charFilename)
606
607
        expIniFilename = self.getProperty("ExpIniFilename").value

608
        self._charTable = ''
609
        if charFilename is None or len(charFilename) <= 0:
610
611
            self.iparmFile = None
            return
612

613
        self._charTable = 'characterizations'
614
615
        results = api.PDLoadCharacterizations(Filename=charFilename,
                                              ExpIniFilename=expIniFilename,
616
                                              OutputWorkspace=self._charTable)
617
        # export the characterizations table
618
        charTable = results[0]
619
620
        if not self.existsProperty("CharacterizationsTable"):
            self.declareProperty(ITableWorkspaceProperty("CharacterizationsTable", self._charTable, Direction.Output))
621
        self.setProperty("CharacterizationsTable", charTable)
622
623

        # get the focus positions from the properties
624
625
626
627
628
629
        self.iparmFile = results[1]
        self._focusPos['PrimaryFlightPath'] = results[2]
        self._focusPos['SpectrumIDs'] = results[3]
        self._focusPos['L2'] = results[4]
        self._focusPos['Polar'] = results[5]
        self._focusPos['Azimuthal'] = results[6]
630

631
    #pylint: disable=too-many-branches
632
    def _load_event_data(self, filename, filter_wall=None, out_ws_name=None, **chunk):
633
634
635
636
637
638
639
640
641
642
643
        """ Load data optionally by chunk strategy
        Purpose:
            Load a complete or partial run, filter bad pulses.
            Output workspace name is formed as
            - user specified (highest priority)
            - instrument-name_run-number_0 (no chunking)
            - instrument-name_run-number_X: X is either ChunkNumber of SpectrumMin
        Requirements:
            1. run number is integer and larger than 0
        Guarantees:
            A workspace is created with name described above
644
        :param filename:
645
        :param filter_wall:
646
647
648
649
650
        :param out_ws_name: name of output workspace specified by user. it will override the automatic name
        :param chunk:
        :return:
        """
        # Check requirements
651
        assert len(filename) > 0, "Input file '%s' does not exist" % filename
652
653
        assert (chunk is None) or isinstance(chunk, dict), 'Input chunk must be either a dictionary or None.'

654
        # get event file's base name
655
656
        base_name = getBasename(filename)

657
658
659
660
661
662
663
664
665
666
667
        # give out the default output workspace name
        if out_ws_name is None:
            if chunk:
                if "ChunkNumber" in chunk:
                    out_ws_name = base_name + "__chk%d" % (int(chunk["ChunkNumber"]))
                elif "SpectrumMin" in chunk:
                    seq_number = 1 + int(chunk["SpectrumMin"])/(int(chunk["SpectrumMax"])-int(chunk["SpectrumMin"]))
                    out_ws_name = base_name + "__chk%d" % seq_number
            else:
                out_ws_name = "%s_%d" % (base_name, 0)
        # END-IF
668

669
        # Specify the other chunk information including Precount, FilterByTimeStart and FilterByTimeStop.
670
        if filename.endswith("_event.nxs") or filename.endswith(".nxs.h5"):
671
            chunk["Precount"] = True
672
673
674
675
676
            if filter_wall is not None:
                if filter_wall[0] > 0.:
                    chunk["FilterByTimeStart"] = filter_wall[0]
                if filter_wall[1] > 0.:
                    chunk["FilterByTimeStop"] = filter_wall[1]
677

678
        # Call Mantid's Load algorithm to load complete or partial data
679
        api.Load(Filename=filename, OutputWorkspace=out_ws_name, **chunk)
680
681

        # Log output
682
        if is_event_workspace(out_ws_name):
683
            self.log().debug("Load run %s: number of events = %d. " % (os.path.split(filename)[-1],
684
                                                                       get_workspace(out_ws_name).getNumberEvents()))
685
        if HAVE_MPI:
686
            msg = "MPI Task = %s ;" % (str(mpiRank))
687
688
            if is_event_workspace(out_ws_name):
                msg += "Number Events = " + str(get_workspace(out_ws_name).getNumberEvents())
689
            self.log().debug(msg)
690

691
        # filter bad pulses
692
        if self._filterBadPulses > 0.:
693
            # record number of events of the original workspace
694
            if is_event_workspace(out_ws_name):
695
                # Event workspace: record original number of events
696
                num_original_events = get_workspace(out_ws_name).getNumberEvents()
697
698
            else:
                num_original_events = -1
699

700
            # filter bad pulse
701
702
            api.FilterBadPulses(InputWorkspace=out_ws_name, OutputWorkspace=out_ws_name,
                                LowerCutoff=self._filterBadPulses)
703

704
            if is_event_workspace(out_ws_name):
705
                # Event workspace
706
                message = "FilterBadPulses reduces number of events from %d to %d (under %s percent) " \
707
708
                          "of workspace %s." % (num_original_events, get_workspace(out_ws_name).getNumberEvents(),
                                                str(self._filterBadPulses), out_ws_name)
709
710
                self.log().information(message)
        # END-IF (filter bad pulse)
711

712
        return out_ws_name
713

714
    def _getStrategy(self, filename):
715
716
        """
        Get chunking strategy by calling mantid algorithm 'DetermineChunking'
717
        :param filename:
718
719
        :return: a list of dictionary.  Each dictionary is a row in table workspace
        """
720
721
        # Determine chunk strategy can search in archive.
        # Therefore this will fail: assert os.path.exists(file_name), 'NeXus file %s does not exist.' % file_name
722
723
724
        self.log().debug("[Fx116] Run file Name : %s,\t\tMax chunk size: %s" %
                         (filename, str(self._chunks)))
        chunks = api.DetermineChunking(Filename=filename, MaxChunkSize=self._chunks)
725

726
        strategy = []
727
728
        for row in chunks:
            strategy.append(row)
729
730
731

        # For table with no rows
        if len(strategy) == 0:
732
            strategy.append({})
733
734
735
736
737

        # delete chunks workspace
        chunks = str(chunks)
        mtd.remove(chunks)

738
739
        return strategy

740
    def __logChunkInfo(self, chunk):
Whitfield, Ross's avatar
Whitfield, Ross committed
741
        keys = sorted(chunk.keys())
742

Zhou, Wenduo's avatar
Zhou, Wenduo committed
743
        keys = [str(key) + "=" + str(chunk[key]) for key in keys]
744
745
        self.log().information("Working on chunk [" + ", ".join(keys) + "]")

746
747
748
    def checkInfoMatch(self, left, right):
        if (left["frequency"].value is not None) and (right["frequency"].value is not None) \
           and (abs(left["frequency"].value - right["frequency"].value)/left["frequency"].value > .05):
Whitfield, Ross's avatar
Whitfield, Ross committed
749
            raise RuntimeError("Cannot add incompatible frequencies (%f!=%f)"
750
751
                               % (left["frequency"].value, right["frequency"].value))
        if (left["wavelength"].value is not None) and (right["wavelength"].value is not None) \
Whitfield, Ross's avatar
Whitfield, Ross committed
752
753
                and abs(left["wavelength"].value - right["wavelength"].value)/left["wavelength"].value > .05:
            raise RuntimeError("Cannot add incompatible wavelengths (%f != %f)"
754
755
                               % (left["wavelength"].value, right["wavelength"].value))

Peterson, Peter's avatar
Peterson, Peter committed
756
    #pylint: disable=too-many-arguments
757
    def _focusAndSum(self, filenames, preserveEvents=True, final_name=None, absorptionWksp=''):
Zhou, Wenduo's avatar
Zhou, Wenduo committed
758
759
760
761
762
763
764
765
766
767
        """Load, sum, and focus data in chunks
        Purpose:
            Load, sum and focus data in chunks;
        Requirements:
            1. input run numbers are in numpy array or list
        Guarantees:
            The experimental runs are focused and summed together
        @param run_number_list:
        @param extension:
        @param preserveEvents:
768
        @param absorptionWksp: will be divided from the data at a per-pixel level
Zhou, Wenduo's avatar
Zhou, Wenduo committed
769
770
        @return: string as the summed workspace's name
        """
771
772
        if final_name is None:
            final_name = getBasename(filenames[0])
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791

        # only pass in the characterizations if the values haven't already been determined
        characterizations = self._charTable
        if self._info is not None:
            characterizations = ''

        # put together a list of the other arguments
        otherArgs = self._focusPos.copy()
        # use the workspaces if they already exists, or pass the filenames down
        # this assumes that AlignAndFocusPowderFromFiles will give the workspaces the canonical names
        cal, grp, msk = [self._instrument + name for name in ['_cal', '_group', '_mask']]
        if mtd.doesExist(cal) and mtd.doesExist(grp) and mtd.doesExist(msk):
            otherArgs['CalibrationWorkspace'] = cal
            otherArgs['GroupingWorkspace'] = grp
            otherArgs['MaskWorkspace'] = msk
        else:
            otherArgs['CalFileName'] = self.calib
            otherArgs['GroupFilename'] = self.getProperty("GroupingFile").value

792
793
        api.AlignAndFocusPowderFromFiles(Filename=','.join(filenames),
                                         OutputWorkspace=final_name,
794
                                         AbsorptionWorkspace=absorptionWksp,
795
796
                                         MaxChunkSize=self._chunks,
                                         FilterBadPulses=self._filterBadPulses,
797
                                         Characterizations=characterizations,
798
                                         CacheDir=self._cache_dir,
799
800
801
802
803
804
                                         Params=self._binning,
                                         ResampleX=self._resampleX,
                                         Dspacing=self._bin_in_dspace,
                                         PreserveEvents=preserveEvents,
                                         RemovePromptPulseWidth=self._removePromptPulseWidth,
                                         CompressTolerance=self.COMPRESS_TOL_TOF,
805
                                         LorentzCorrection=self._lorentz,
806
807
808
809
810
                                         UnwrapRef=self._LRef,
                                         LowResRef=self._DIFCref,
                                         LowResSpectrumOffset=self._lowResTOFoffset,
                                         CropWavelengthMin=self._wavelengthMin,
                                         CropWavelengthMax=self._wavelengthMax,
811
812
                                         FrequencyLogNames=self.getProperty("FrequencyLogNames").value,
                                         WaveLengthLogNames=self.getProperty("WaveLengthLogNames").value,