PeaksOnSurface

[Owen Arnold]

This is an adapted version of PeaksInRegion. Which is better suitable for the SliceViewer peaks overlay.

• We only actually need to determine Peaks interaction with 1 plane face, rather than 6 because the slice plane is infinitely thin. This should cut down a large amount of wasted computation.
• The CheckPeakExtents option is always on in the SliceViewer, this is consistent with how we would want to work with a plane.
• We can return the absolute distance between the peak center and the plane. This is information we currently throw away in PeaksInRegion, and then recalculate separately in the SliceViewer. Wouldn't it just be better to have an additional column in the output peaks workspace to contain this information.

We should leave PeaksInRegion in place, because this could be quite a handy algorithm for some future work.

I'm going to wire this algorithm in under a separate ticket, to avoid this ticket bloating out. Also, because it's going to be very hard to validate both the algorithm itself, and it's new usage/consumption in one-shot.

[Owen Arnold]

Tester: Run this

Load(Filename='SXD23767.raw',OutputWorkspace='SXD23767',LoadMonitors='Exclude')
ConvertToDiffractionMDWorkspace(InputWorkspace='SXD23767',OutputWorkspace='QLab',LorentzCorrection='1',SplitThreshold='50',MaxRecursionDepth='13',Extents='-15,15,-15,15,-15,15')
peaks_qLab = FindPeaksMD(InputWorkspace='QLab',PeakDistanceThreshold='1',MaxPeaks='60')
peakRadius = 0.1;
detectorZPlanePosition =  -0.278  # For bank 11, all detectors sit on constant y.
yPlane = detectorZPlanePosition - peakRadius
queryPeaksInRegion = PeaksInRegion(CheckPeakExtents=True,InputWorkspace=peaks_qLab,CoordinateFrame='Detector space',PeakRadius=peakRadius,Extents=[-0.09,0.09,yPlane,yPlane,-0.09,0.09])
queryPeaksOnSurface = PeaksOnSurface(InputWorkspace=peaks_qLab,CoordinateFrame='Detector space',PeakRadius=peakRadius,Vertex1=[-0.09,yPlane,-0.09],Vertex2=[-0.09,yPlane,0.09],Vertex3=[0.09,yPlane,0.09],Vertex4=[0.09,yPlane,-0.09])

The outputs queryPeakInRegion and queryPeaksOnSurface should be the same (ignore the distance column output).

Use the script above as the basis for running other comparison tests, for example, try increasing and reducing the radius parameter.

The algorithm PeaksOnSurface will also work in a non-axis aligned fashion. It would be a good idea to test this out too.

Please describe what tests you have run in the comments of this ticket.

[Owen Arnold]

refs #7325. Move working into base class.

New algorithm base class PeaksIntersection defined for algorithms like PeaksInRegion and any other new algorithm with surfaces involved. Had to make some methods virtual in the process of moving working from PeaksInRegion into the base class, but this does not appear to have affected performance. Will still process 40000 peaks in about 0.3 of a second.

Changeset: aa951b2a5c0c85e0b19414f21a6f0cd80de5801b

[Owen Arnold]

refs #7325. Implement PeaksOnSurface algorithm.

Changeset: 9fc7daf6ea9e4261c4040e232bdd3072ba3fba0d

[Owen Arnold]

refs #7325. Add performance tests.

Comparing against peaks in region, the equivalent performance tests are producing similar results. PeaksOnSurface is very slightly slower. But we can improve performance here yet. Timings are 0.039 for PeaksInRegion and 0.042 for PeaksOnSurface.

Changeset: f5f331a714e182cfa686f48242ec53c300af1cc1

[Owen Arnold]

refs #7325. Document geometric calculations.

Changeset: d5de466d6a3b7a76d5d49f44c42ab902b8b905ad

[Owen Arnold]

refs #7325. Extend output workspace.

Changeset: 2db16c59decafc9d394198c303bfb87631e4d96a

[Owen Arnold]

refs #7325. Remove property from PeaksOnSurface.

Changeset: 29aa2da9b8c6b301826f707de9f9ed1fca552ea5

[Owen Arnold]

refs #7325. Explicitly specify container type.

Changeset: 9c9569502efdc534af78e26e1c074530e471ef76

[Owen Arnold]

refs #7325. Fix type conversion warning.

Changeset: 139282f6930976367ca1a48d0803d773291f7223

[Owen Arnold]

branch is feature/7325_peaks_on_surface

[Owen Arnold]

I need to investigate the following script provided by the tester:

run_segfault = False

peakRadius = 0.01;
detectorZPlanePosition =  -0.278  # For bank 11, all detectors sit on constant y.
yPlane = detectorZPlanePosition - peakRadius

Load(Filename='SXD23767.raw', OutputWorkspace='SXD23767',  LoadMonitors='Exclude')

if not run_segfault:
    ConvertToDiffractionMDWorkspace(InputWorkspace='SXD23767', OutputWorkspace='QLab', LorentzCorrection='1', SplitThreshold='50', MaxRecursionDepth='13', Extents='-15,15,-15,15,-15,15')

    # Trying to reorient the axes via basis vectors. This works but only finds 4 
    # peaks instead of 60. Also, FindPeaksMD won't work with non-standard axis
    # labels, so I had to co-op the standard ones.

    BinMD(InputWorkspace='QLab', AxisAligned='0', BasisVector0='Q_lab_x,A,1,1,0', BasisVector1='Q_lab_y,A,-1,1,0', BasisVector2='Q_lab_z,A,0,0,1', OutputExtents='-15,15,-15,15,-15,15', OutputBins='10,10,10', Parallel='1', OutputWorkspace='QLab_Rotb')

    peaks_QLab_Rotb = FindPeaksMD(InputWorkspace='QLab_Rotb', PeakDistanceThreshold='1', MaxPeaks='60')

    queryPeaksInRegion_b = PeaksInRegion(CheckPeakExtents=True, InputWorkspace=peaks_QLab_Rotb, CoordinateFrame='Detector space',   PeakRadius=peakRadius, Extents=[-0.09,0.09,yPlane,yPlane,-0.09,0.09])

    queryPeaksOnSurface_b = PeaksOnSurface(InputWorkspace=peaks_QLab_Rotb, CoordinateFrame='Detector space', PeakRadius=peakRadius, Vertex1=[-0.09,yPlane,-0.09], Vertex2=[-0.09,yPlane,0.09], Vertex3=[0.09,yPlane,0.09], Vertex4=[0.09,yPlane,-0.09])

else:
    # Trying to do reorientation at convert time. This finds 60 peaks, but causes 
    # segfault in PeaksInRegion.

    ConvertToMD(InputWorkspace='SXD23767', OutputWorkspace='QLab_Rotc', QDimensions='Q3D', dEAnalysisMode='Elastic', Q3DFrames='Q_lab', LorentzCorrection='1', MinValues='-15,-15,-15', MaxValues='15,15,15', Uproj='1,1,0', Vproj='-1,1,0', Wproj='0,0,1', SplitInto='2', SplitThreshold='50', MaxRecursionDepth='13', MinRecursionDepth='5')

    peaks_QLab_Rotc = FindPeaksMD(InputWorkspace='QLab_Rotc', PeakDistanceThreshold='1', MaxPeaks='60')

    queryPeaksInRegion_c = PeaksInRegion(CheckPeakExtents=True, InputWorkspace=peaks_QLab_Rotc, CoordinateFrame='Detector space', PeakRadius=peakRadius, Extents=[-0.09,0.09,yPlane,yPlane,-0.09,0.09])

    queryPeaksOnSurface_c = PeaksOnSurface(InputWorkspace=peaks_QLab_Rotc, CoordinateFrame='Detector space', PeakRadius=peakRadius, Vertex1=[-0.09,yPlane,-0.09], Vertex2=[-0.09,yPlane,0.09], Vertex3=[0.09,yPlane,0.09], Vertex4=[0.09,yPlane,-0.09])

[Owen Arnold]

OK, so the second part of this is a genuine crash owing to the fact that it's trying to access the location of a detector that is null. This does need to be fixed as part of this ticket.

I'm confused about the first part (of the script above) about which part appears to be wrong. Slightly modifying the script (see below). The number of peaks found after rotation is only 4, while 60 are found before rotation. Also, the act of binning the data to a regular grid, will ruin the ability of FindPeaksMD to find peaks, given that it's the natural adaptive binning structure, which is particularly useful for SC peak finding.

Not sure how to proceed from here.

run_segfault = False

peakRadius = 0.01;
detectorZPlanePosition =  -0.278  # For bank 11, all detectors sit on constant y.
yPlane = detectorZPlanePosition - peakRadius

Load(Filename='SXD23767.raw', OutputWorkspace='SXD23767',  LoadMonitors='Exclude')

if not run_segfault:
    ConvertToDiffractionMDWorkspace(InputWorkspace='SXD23767', OutputWorkspace='QLab', LorentzCorrection='1', SplitThreshold='50', MaxRecursionDepth='13', Extents='-15,15,-15,15,-15,15')

    # Trying to reorient the axes via basis vectors. This works but only finds 4 
    # peaks instead of 60. Also, FindPeaksMD won't work with non-standard axis
    # labels, so I had to co-op the standard ones.
    
    peaks_QLab_Rotb_prior = FindPeaksMD(InputWorkspace='QLab', PeakDistanceThreshold='1', MaxPeaks='60')

    BinMD(InputWorkspace='QLab', AxisAligned='0', BasisVector0='Q_lab_x,A,1,1,0', BasisVector1='Q_lab_y,A,-1,1,0', BasisVector2='Q_lab_z,A,0,0,1', OutputExtents='-15,15,-15,15,-15,15', OutputBins='10,10,10', Parallel='1', OutputWorkspace='QLab_Rotb')

    peaks_QLab_Rotb = FindPeaksMD(InputWorkspace='QLab_Rotb', PeakDistanceThreshold='1', MaxPeaks='60')

[Michael Reuter]

I would say that if the first case is detrimental to SCD data, then no one will be using that procedure. Therefore, I'm fine with letting that one go and not considering it any further.

[Owen Arnold]

Now that #7386 is done, re-execute the script with run_segfault=True, and the algorithm should fall over, but not crash MantidPlot. This is sufficient because the workflow to get to the point of this algorithm failing is actually unrealistic for any kind of single crystal diffraction experiments.

[Michael Reuter]

OK, the script does not crash Mantid anymore and the original script works as before. For the segfault case, the resulting error message

PeaksInRegion: Detector cannot be found. search object Peak at row -1
Error in execution of algorithm PeaksInRegion:
PeaksInRegion: error (see log)

especially the last line, leads one to believe there is more information to be had. I can't find any more information anywhere. I think the error message needs to be clarified a bit to either add the extra information or remove the extraneous message. I've opened #7400 to address this.

[Stuart Campbell]

This ticket has been transferred to github issue ​8171