Ticket #6: calculate_monitor_M3toM2_and_IQ_empty.py

# reduction without sensitivity and efficiency
from mantid.simpleapi import *
import math

def mergeWS(prepath, wsList,output):
    for ix in wsList:
        name='RUN'+str(ix).zfill(7)
        LoadNexus(Filename=prepath+'\\'+name+'\\'+output+'.nxs',OutputWorkspace=str(ix))
    for i in range(len(wsList)):
        if len(wsList)== 1:
            CloneWorkspace(InputWorkspace=str(wsList[i]), OutputWorkspace=output)
            DeleteWorkspace(Workspace=str(wsList[i]))
        elif wsList[i+1]!=wsList[-1]:
            Plus(LHSWorkspace=str(wsList[i]),RHSWorkspace=str(wsList[i+1]), OutputWorkspace=str(wsList[i+1]))
            DeleteWorkspace(Workspace=str(wsList[i]))
        else:
            Plus(LHSWorkspace=str(wsList[i]),RHSWorkspace=str(wsList[i+1]), OutputWorkspace=output)
            DeleteWorkspace(Workspace=str(wsList[i]))
            DeleteWorkspace(Workspace=str(wsList[i+1]))
            break


def getMask(xPids, yPids):
    index=[]
    xnum=xPids[1]-xPids[0]+1
    ynum=yPids[1]-yPids[0]+1
    sta=yPids[0]*125+xPids[0]
    for i in range(xnum):
        for j in range(ynum):
            index.append(sta+j*125+i)
    print index
    return index


def findCmaskWS(wsname):
    x=0.0009
    y=-0.00229
    MoveInstrumentComponent(Workspace=wsname, ComponentName='Bank01', X=x, Y=y, RelativePosition=False)
    index=getMask([0,18],[0,119])
    MaskDetectors(Workspace=wsname, WorkspaceIndexList=index)
    index=getMask([105,124],[0,119])
    MaskDetectors(Workspace=wsname, WorkspaceIndexList=index)
    inner=0.048
    xml_inner='<infinite-cylinder id="beam_stop"><centre x="0" y="0" z="0.0" /><axis x="0" y="0" z="1" /><radius val= "'+str(inner)+'" /></infinite-cylinder><algebra val="beam_stop"/>'
    MaskDetectorsInShape(Workspace=wsname, ShapeXML=xml_inner)

    index=getMask([60,63],[0,57])
    MaskDetectors(Workspace=wsname, WorkspaceIndexList=index)
    
    
    

empty_run=[960,961]

monSUM=True
prepath=r'C:\Users\admin\Documents\python\CSNSSANS'


#empty
mergeWS(prepath, empty_run, 'sample')
CloneWorkspace(InputWorkspace='sample',OutputWorkspace='empty_tmp')

LoadInstrument(Workspace='empty_tmp', Filename=prepath+'/paramData/detector_mod.xml', RewriteSpectraMap=True)
#CloneWorkspace(InputWorkspace='bg',OutputWorkspace='empty_tmp')
findCmaskWS('empty_tmp')
ConvertUnits(InputWorkspace='empty_tmp', OutputWorkspace='empty_tmp0', Target='Wavelength', AlignBins=True)


wave_rebin='2,0.02,7.9'

#load calibration 
#==========================================
#LoadAscii(Filename=prepath+'/paramData/cal_1M.dat', OutputWorkspace='cal', Separator='Automatic', Unit='Wavelength')
LoadAscii(Filename=prepath+'/paramData/cal_bp_1M.dat', OutputWorkspace='cali', Separator='Automatic', Unit='Wavelength')
#ConvertToHistogram(InputWorkspace='cal', OutputWorkspace='cal')

#LoadNexus(Filename=prepath+'\paramData\cal_1M.nxs', OutputWorkspace='cal')
ConvertToHistogram(InputWorkspace='cali', OutputWorkspace='cali')
Rebin(InputWorkspace='cali', OutputWorkspace='cali', Params=wave_rebin)
#Rebin(InputWorkspace='cal', OutputWorkspace='cal', Params=wave_rebin)
#Divide(LHSWorkspace="cal", RHSWorkspace='v', OutputWorkspace="cal", AllowDifferentNumberSpectra=True)


#empty
'''
mergeWS(prepath, empty_run, 'sample')
CloneWorkspace(InputWorkspace='sample',OutputWorkspace='empty_tmp')
LoadInstrument(Workspace='empty_tmp', Filename=prepath+'/paramData/detector_mod.xml', RewriteSpectraMap=True)
#CloneWorkspace(InputWorkspace='bg',OutputWorkspace='empty_tmp')
findCmaskWS('empty_tmp')
ConvertUnits(InputWorkspace='empty_tmp', OutputWorkspace='empty_tmp', Target='Wavelength', AlignBins=True)

'''
Rebin(InputWorkspace='empty_tmp0',OutputWorkspace='empty_tmp',Params=wave_rebin)

mergeWS(prepath, empty_run,'monitor2')
RenameWorkspace(InputWorkspace='monitor2',OutputWorkspace='mon_empty')
#CropWorkspace(InputWorkspace="mon_empty", OutputWorkspace="empty_m2",StartWorkspaceIndex=0, EndWorkspaceIndex=0)
ConvertUnits(InputWorkspace='mon_empty', OutputWorkspace='empty_m2', Target='Wavelength', AlignBins=True)
Rebin(InputWorkspace='empty_m2',OutputWorkspace='empty_m2',Params=wave_rebin)
#RebinToWorkspace(WorkspaceToRebin='empty_m2',WorkspaceToMatch='empty_tmp',OutputWorkspace='empty_m2')

#Divide(LHSWorkspace="empty_m2", RHSWorkspace='cal', OutputWorkspace="empty_m2", AllowDifferentNumberSpectra=True)
#'''
if monSUM:
    name=mtd['empty_m2']
    counts=sum(name.readY(0))*1.0
    CreateSingleValuedWorkspace(OutputWorkspace='m2',DataValue=counts)
    Divide(LHSWorkspace='empty_tmp', RHSWorkspace='m2', OutputWorkspace='empty_tmp', AllowDifferentNumberSpectra=True)
else:
    Divide(LHSWorkspace='empty_tmp', RHSWorkspace='empty_m2', OutputWorkspace='empty_tmp', AllowDifferentNumberSpectra=True)

Divide(LHSWorkspace='empty_tmp', RHSWorkspace='cali', OutputWorkspace='empty_tmp', AllowDifferentNumberSpectra=True)
#Multiply(LHSWorkspace='empty_tmp', RHSWorkspace='cali', OutputWorkspace='empty_tmp', AllowDifferentNumberSpectra=True)


q_rebin='0.01,-0.08,0.4'
#q_rebin='0.01,0.01,0.4'
#q_rebin='0.01,0.005,0.8'
Q1D(DetBankWorkspace='empty_tmp', OutputWorkspace='IQ_air', OutputBinning=q_rebin, AccountForGravity=False, SolidAngleWeighting=True)
ReplaceSpecialValues(InputWorkspace='IQ_air', OutputWorkspace='IQ_air', NaNValue=0, InfinityValue=0)


dlamda=0.5
CreateSingleValuedWorkspace(OutputWorkspace='dlamda_sample', DataValue=dlamda)
Multiply(LHSWorkspace='IQ_air', RHSWorkspace='dlamda_sample', OutputWorkspace='IQ_air', AllowDifferentNumberSpectra=True)


trans_run=[1118]
mergeWS(prepath, trans_run,'monitors')
RenameWorkspace(InputWorkspace='monitors',OutputWorkspace='monTran')




LoadInstrument(Workspace='monTran', Filename=prepath+'/paramData/monitors.xml', RewriteSpectraMap=True)
ConvertUnits(InputWorkspace='monTran', OutputWorkspace='monTran', Target='Wavelength', AlignBins=True)


Rebin(InputWorkspace='monTran', OutputWorkspace="monTran", Params=wave_rebin)

name=mtd['monTran']
monTranM2=(name.readY(0))*1.0
monTranM3=(name.readY(1))*1.0
monTranMWavelength=(name.readX(0))*1.0

CreateWorkspace(OutputWorkspace='monTranM3',DataX=monTranMWavelength,DataY=monTranM3)
CreateWorkspace(OutputWorkspace='monTranM2',DataX=monTranMWavelength,DataY=monTranM2)

Divide(LHSWorkspace='monTranM3', RHSWorkspace='monTranM2',OutputWorkspace='monTranM3toM2'+str(trans_run), AllowDifferentNumberSpectra=True)







































', OutputWorkspace='empty_tmp', AllowDifferentNumberSpectra=True)