Context Navigation

Back to Ticket #7670

Ticket #7670: _procedures.py

File _procedures.py, 59.8 KB (added by Gesner Passos, 7 years ago)
The file for the new api.

Line
1	from math import *
2	from mantid.simpleapi import *
3	from mantid.api import WorkspaceGroup
4	try:
5	from mantidplot import *
6	except ImportError:
7	pass
8
9	import numpy as n
10
11	def addRuns(runlist,wname):
12	#DeleteWorkspace(str(wname))
13	output=str(wname)
14	if runlist[0] != "0":
15	#nzeros=8-len(str(runlist[0]))
16	#fpad=""
17	#for i in range(nzeros):
18	# fpad+="0"
19	#filename="offspec"+fpad+str(runlist[0])+".nxs"
20	#fname=str(FileFinder.findRuns(filename))
21	#fname=str.replace(fname,'\\','/')
22	#fname=str.replace(fname,'[','')
23	#fname=str.replace(fname,']','')
24	#fname=str.replace(fname,'\'','')
25	#fname=fname.lower()
26	##fname=str.replace(fname,'.nxs','.raw')
27	#Load(fname,output)
28	Load(str(runlist[0]),OutputWorkspace=output)
29	else:
30	#dae="ndx"+config['default.instrument'].lower()
31	dae="ndxoffspec"
32	LoadDAE(DAEname=dae,OutputWorkspace=output,SpectrumMin="1")
33	#LoadLiveData(Instrument="OFFSPEC",AccumulationMethod="Replace",OutputWorkspace="output")
34	if isinstance(mtd[output], WorkspaceGroup):
35	for k in mtd[output].getNames():
36	mtd[k].setYUnit('Counts')
37	else:
38	mtd[output].setYUnit('Counts')
39
40	if len(runlist) > 1:
41	for i in range(1,len(runlist)):
42	if runlist[i] != "0":
43	#nzeros=8-len(str(runlist[i]))
44	#fpad=""
45	#for j in range(nzeros):
46	# fpad+="0"
47	#filename="offspec"+fpad+str(runlist[i])+".nxs"
48	#fname=str(FileFinder.findRuns(filename))
49	#fname=str.replace(fname,'\\','/')
50	#fname=str.replace(fname,'[','')
51	#fname=str.replace(fname,']','')
52	#fname=str.replace(fname,'\'','')
53	#fname=fname.lower()
54	##fname=str.replace(fname,'.nxs','.raw')
55	#Load(fname,"wtemp")
56	Load(str(runlist[i]),OutputWorkspace="wtemp")
57	else:
58	#dae="ndx"+config['default.instrument'].lower()
59	dae="ndxoffspec"
60	LoadDAE(DAEname=dae,OutputWorkspace="wtemp",SpectrumMin="1")
61	#LoadLiveData(Instrument="OFFSPEC",AccumulationMethod="Replace",OutputWorkspace="output")
62	if isinstance(mtd['wtemp'], WorkspaceGroup):
63	for k in mtd['wtemp'].getNames():
64	mtd[k].setYUnit('Counts')
65	else:
66	mtd[output].setYUnit('Counts')
67	Plus(output,"wtemp",OutputWorkspace=output)
68	DeleteWorkspace("wtemp")
69
70	#logger.notice("addRuns Completed")
71	#
72	#===================================================================================================================
73	#
74	'''
75	parse a text string of the format "1-6:2+8+9,10+11+12+13-19:3,20-24"
76	to return a structure containing the separated lists [1, 3, 5, 8, 9],
77	[10, 11, 12, 13, 16, 19] and [20, 21, 22, 23, 24]
78	as integer lists that addRuns can handle.
79	'''
80	def parseRunList(istring):
81	if len(istring) >0:
82	s1=istring.split(',')
83	rlist1=[]
84	for i in range(len(s1)):
85	tstr=s1[i].strip()
86	if len(tstr) > 0:
87	rlist1.append(tstr)
88	rlist=[]
89	for i in range(len(rlist1)):
90	rlist2=[]
91	if rlist1[i].find('+') >= 0:
92	tstr=rlist1[i].split('+')
93	for j in range(len(tstr)):
94	if tstr[j].find(':') >=0 and tstr[j].find('-') >=0:
95	tstr[j].strip()
96	tstr2=tstr[j].split('-')
97	tstr3=tstr2[1].split(':')
98	r1=range(int(tstr2[0]),int(tstr3[0])+1,int(tstr3[1]))
99	for k in r1:
100	rlist2.append(str(k))
101	elif tstr[j].find('-') >=0:
102	tstr[j].strip()
103	tstr2=tstr[j].split('-')
104	r1=range(int(tstr2[0]),int(tstr2[1])+1)
105	for k in r1:
106	rlist2.append(str(k))
107	else:
108	rlist2.append(tstr[j])
109	else:
110	if rlist1[i].find(':') >=0 and rlist1[i].find('-')>=0:
111	rlist1[i].strip()
112	tstr2=rlist1[i].split('-')
113	tstr3=tstr2[1].split(':')
114	r1=range(int(tstr2[0]),int(tstr3[0])+1,int(tstr3[1]))
115	for k in r1:
116	rlist2.append(str(k))
117	elif rlist1[i].find('-')>=0:
118	rlist1[i].strip()
119	tstr2=rlist1[i].split('-')
120	r1=range(int(tstr2[0]),int(tstr2[1])+1)
121	for k in r1:
122	rlist2.append(str(k))
123	else:
124	rlist2.append(rlist1[i])
125	rlist.append(rlist2)
126	return rlist
127	#
128	#===================================================================================================================
129	#
130	def parseNameList(istring):
131	s1=istring.split(',')
132	namelist=[]
133	for i in range(len(s1)):
134	tstr=s1[i].strip()
135	namelist.append(tstr)
136	return namelist
137
138	#
139	#===================================================================================================================
140	#
141	def floodnorm(wkspName,floodfile):
142	#
143	# pixel by pixel efficiency correction for the linear detector
144	#
145	flood_wksp = "ld240flood"
146	if flood_wksp not in mtd:
147	LoadNexusProcessed(Filename=floodfile,OutputWorkspace=flood_wksp)
148
149	Divide(LHSWorkspace=wkspName, RHSWorkspace=flood_wksp, OutputWorkspace=wkspName)
150	#
151	#===================================================================================================================
152	#
153	# Plot a bunch of workspaces as 2D maps
154	# using the supplied limits and log scale settings
155	#
156	def plot2D(wkspNames,limits,logScales):
157	nplot=0
158	workspace_mtx=[]
159	wNames=parseNameList(wkspNames)
160	for i in range(len(wNames)):
161	w1=mtd[wNames[i]]
162	if isinstance(w1, WorkspaceGroup):
163	w1names=w1.getNames()
164	for j in range(len(w1names)):
165	#workspace_mtx.append(mantidplot.importMatrixWorkspace(w1names[j]))
166	workspace_mtx.append(mantidplot.importMatrixWorkspace(w1names[j]))
167	gr2d=workspace_mtx[nplot].plotGraph2D()
168	nplot=nplot+1
169	l=gr2d.activeLayer()
170	if logScales[0]=="0":
171	l.setAxisScale(mantidplot.Layer.Bottom,float(limits[0]),float(limits[1]))
172	elif logScales[0]=="2":
173	l.setAxisScale(mantidplot.Layer.Bottom,float(limits[0]),float(limits[1]),1)
174	if logScales[1]=="0":
175	l.setAxisScale(mantidplot.Layer.Left,float(limits[2]),float(limits[3]))
176	elif logScales[1]=="2":
177	l.setAxisScale(mantidplot.Layer.Left,float(limits[2]),float(limits[3]),1)
178	if logScales[2]=="0":
179	l.setAxisScale(mantidplot.Layer.Right,float(limits[4]),float(limits[5]))
180	elif logScales[2]=="2":
181	l.setAxisScale(mantidplot.Layer.Right,float(limits[4]),float(limits[5]),1)
182	else:
183	workspace_mtx.append(mantidplot.importMatrixWorkspace(wNames[i]))
184	gr2d=workspace_mtx[nplot].plotGraph2D()
185	nplot=nplot+1
186	l=gr2d.activeLayer()
187	if logScales[0]=="0":
188	l.setAxisScale(mantidplot.Layer.Bottom,float(limits[0]),float(limits[1]))
189	elif logScales[0]=="2":
190	l.setAxisScale(mantidplot.Layer.Bottom,float(limits[0]),float(limits[1]),1)
191	if logScales[1]=="0":
192	l.setAxisScale(mantidplot.Layer.Left,float(limits[2]),float(limits[3]))
193	elif logScales[1]=="2":
194	l.setAxisScale(mantidplot.Layer.Left,float(limits[2]),float(limits[3]),1)
195	if logScales[2]=="0":
196	l.setAxisScale(mantidplot.Layer.Right,float(limits[4]),float(limits[5]))
197	elif logScales[2]=="2":
198	l.setAxisScale(mantidplot.Layer.Right,float(limits[4]),float(limits[5]),1)
199
200	logger.notice("plot2D finished")
201	#
202	#===================================================================================================================
203	#
204	# Plot a bunch of workspaces as 2D maps
205	# using the supplied limits and log scale settings
206	#
207	def XYPlot(wkspNames,spectra,limits,logScales,errors,singleFigure):
208	wNames=parseNameList(wkspNames)
209	spec=parseNameList(spectra)
210	ploterr=0
211	xLog=0
212	yLog=0
213	if errors == "2":
214	ploterr=1
215	if logScales[0] == "2":
216	xLog=1
217	if logScales[1] == "2":
218	yLog=1
219
220	if singleFigure == "2":
221	p1=plotSpectrum(wNames,spec,ploterr)
222	l=p1.activeLayer()
223	l.setAxisScale(mantidplot.Layer.Bottom,float(limits[0]),float(limits[1]),xLog)
224	l.setAxisScale(mantidplot.Layer.Left,float(limits[2]),float(limits[3]),yLog)
225	else:
226	for i in range(len(wNames)):
227	p1=plotSpectrum(wNames[i],spec,ploterr)
228	l=p1.activeLayer()
229	l.setAxisScale(mantidplot.Layer.Bottom,float(limits[0]),float(limits[1]),xLog)
230	l.setAxisScale(mantidplot.Layer.Left,float(limits[2]),float(limits[3]),yLog)
231
232	logger.notice("XYPlot finished")
233	#
234	#===================================================================================================================
235	#
236	def nrtestfn(runlist,wnames):
237
238	rlist=parseRunList(runlist)
239	logger.notice("This is the runlist:"+str(rlist))
240	namelist=parseNameList(wnames)
241	logger.notice("This is the nameslist:"+str(namelist))
242	for i in range(len(rlist)):
243	addRuns(rlist[i],namelist[i])
244	#Load(Filename="L:/RawData/cycle_10_1/OFFSPEC0000"+str(rlist[i][j])+".nxs",OutputWorkspace="w"+str(rlist[i][j]),LoaderName="LoadISISNexus")
245
246	logger.notice("nrtestfn completed")
247	'''
248	output="w7503"
249	plotper=[1,2]
250	rebpars="0.5,0.025,14.5"
251	spmin=0
252	spmax=239
253	Load(Filename="L:/RawData/cycle_10_1/OFFSPEC00007503.nxs",OutputWorkspace=output,LoaderName="LoadISISNexus")
254	workspace_mtx=[]
255	nplot=0
256	for i in plotper:
257	ConvertUnits(InputWorkspace=output+"_"+str(i),OutputWorkspace=output+"_"+str(i),Target="Wavelength",AlignBins="1")
258	Rebin(InputWorkspace=output+"_"+str(i),OutputWorkspace=output+"_"+str(i),Params=rebpars)
259	CropWorkspace(InputWorkspace=output+"_"+str(i),OutputWorkspace=output+"_"+str(i)+"m",StartWorkspaceIndex="1",EndWorkspaceIndex="1")
260	CropWorkspace(InputWorkspace=output+"_"+str(i),OutputWorkspace=output+"_"+str(i)+"d",StartWorkspaceIndex=str(spmin+4),EndWorkspaceIndex=str(spmax+4))
261	Divide(output+"_"+str(i)+"d",output+"_"+str(i)+"m",OutputWorkspace=output+"_"+str(i)+"n")
262	workspace_mtx.append(mantidplot.importMatrixWorkspace(output+"_"+str(i)+"n"))
263	gr2d=workspace_mtx[nplot].plotGraph2D()
264	nplot=nplot+1
265	l=gr2d.activeLayer()
266
267	logger.notice("quickPlot Finished")
268	'''
269	def removeoutlayer(wksp):
270	# remove counts from bins where there are so few counts it makes a mess of the polarisation
271	# calculation
272	a1=mtd[wksp]
273	nspec=a1.getNumberHistograms()
274	x=a1.readX(0)
275	for i in range(nspec):
276	for j in range(len(x)-1):
277	y=a1.readY(i)[j]
278	if (y<2):
279	a1.dataY(i)[j]=0.0;
280	a1.dataE(i)[j]=0.0;
281
282	def nrSESANSFn(runList,nameList,P0runList,P0nameList,minSpec,maxSpec,upPeriod,downPeriod,existingP0,SEConstants,gparams,convertToSEL,lnPOverLam,diagnostics="0",removeoutlayer="0",floodfile="none",):
283	nlist=parseNameList(nameList)
284	logger.notice("This is the sample nameslist:"+str(nlist))
285	rlist=parseRunList(runList)
286	logger.notice("This is the sample runlist:"+str(rlist))
287	for i in range(len(rlist)):
288	addRuns(rlist[i],nlist[i])
289
290	P0nlist=parseNameList(P0nameList)
291	logger.notice("This is the P0nameslist:"+str(P0nlist))
292	if existingP0 != "2":
293	P0rlist=parseRunList(P0runList)
294	logger.notice("This is the P0runlist:"+str(P0rlist))
295	for i in range(len(P0rlist)):
296	addRuns(P0rlist[i],P0nlist[i])
297
298	mon_spec=int(gparams[3])-1
299	minSp=int(minSpec)-1
300	maxSp=int(maxSpec)-1
301	reb=gparams[0]+","+gparams[1]+","+gparams[2]
302	if len(gparams) == 5:
303	mapfile=gparams[4]
304
305	for i in nlist:
306	a1=mtd[i+"_1"]
307	nspec=a1.getNumberHistograms()
308	if nspec == 1030 and minSp != 3:
309	GroupDetectors(InputWorkspace=i,OutputWorkspace=i,MapFile=mapfile)
310	ConvertUnits(InputWorkspace=i,OutputWorkspace=i,Target="Wavelength",AlignBins=1)
311	Rebin(i,reb,OutputWorkspace=i)
312	if (removeoutlayer != "0"):
313	removeoutlayer(i+"_1")
314	removeoutlayer(i+"_2")
315	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
316	if nspec == 245:
317	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"2ddet",StartWorkspaceIndex=4,EndWorkspaceIndex=243)
318	if (floodfile != "none"):
319	floodnorm(i+"2ddet",floodfile)
320	if nspec == 1030:
321	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"2ddet",StartWorkspaceIndex=3,EndWorkspaceIndex=124)
322	if int(maxSpec) > int(minSpec):
323	SumSpectra(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=minSp,EndWorkspaceIndex=maxSp)
324	else:
325	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=minSp,EndWorkspaceIndex=maxSp)
326
327	Divide(LHSWorkspace=i+"det",RHSWorkspace=i+"mon",OutputWorkspace=i+"norm")
328	if nspec > 4 and minSp != 3:
329	Divide(LHSWorkspace=i+"2ddet",RHSWorkspace=i+"mon",OutputWorkspace=i+"2dnorm")
330	DeleteWorkspace(i+"mon")
331	if (diagnostics == "0"):
332	DeleteWorkspace(i+"det")
333	DeleteWorkspace(i)
334	Minus(LHSWorkspace=i+"norm_"+upPeriod,RHSWorkspace=i+"norm_"+downPeriod,OutputWorkspace="num")
335	Plus(LHSWorkspace=i+"norm_2",RHSWorkspace=i+"norm_1",OutputWorkspace="den")
336	Divide(LHSWorkspace="num",RHSWorkspace="den",OutputWorkspace=i+"pol")
337	ReplaceSpecialValues(InputWorkspace=i+"pol",OutputWorkspace=i+"pol",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
338	if nspec >4 and minSp != 3:
339	#print i+"2dnorm_"+upPeriod
340	#print i+"2dnorm_"+downPeriod
341	Minus(LHSWorkspace=i+"2dnorm_"+upPeriod,RHSWorkspace=i+"2dnorm_"+downPeriod,OutputWorkspace="num")
342	Plus(LHSWorkspace=i+"2dnorm_2",RHSWorkspace=i+"2dnorm_1",OutputWorkspace="den")
343	Divide(LHSWorkspace="num",RHSWorkspace="den",OutputWorkspace=i+"2dpol")
344	ReplaceSpecialValues(InputWorkspace=i+"2dpol",OutputWorkspace=i+"2dpol",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
345	#DeleteWorkspace(i+"norm_2")
346	#DeleteWorkspace(i+"norm_1")
347	DeleteWorkspace("num")
348	DeleteWorkspace("den")
349
350	if existingP0 != "2":
351	for i in P0nlist:
352	ConvertUnits(InputWorkspace=i,OutputWorkspace=i,Target="Wavelength",AlignBins=1)
353	Rebin(InputWorkspace=i,OutputWorkspace=i,Params=reb)
354	removeoutlayer(i+"_1")
355	removeoutlayer(i+"_2")
356	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
357	if int(maxSpec) > int(minSpec):
358	SumSpectra(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=minSp,EndWorkspaceIndex=maxSp)
359	else:
360	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=minSp,EndWorkspaceIndex=maxSp)
361	Divide(LHSWorkspace=i+"det",RHSWorkspace=i+"mon",OutputWorkspace=i+"norm")
362	DeleteWorkspace(i+"mon")
363	DeleteWorkspace(i+"det")
364	DeleteWorkspace(i)
365	Minus(LHSWorkspace=i+"norm_"+upPeriod,RHSWorkspace=i+"norm_"+downPeriod,OutputWorkspace="num")
366	Plus(LHSWorkspace=i+"norm_2",RHSWorkspace=i+"norm_1",OutputWorkspace="den")
367	Divide(LHSWorkspace="num",RHSWorkspace="den",OutputWorkspace=i+"pol")
368	ReplaceSpecialValues(InputWorkspace=i+"pol",OutputWorkspace=i+"pol",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
369	DeleteWorkspace(i+"norm_2")
370	DeleteWorkspace(i+"norm_1")
371	DeleteWorkspace("num")
372	DeleteWorkspace("den")
373
374	for i in range(len(nlist)):
375	if existingP0 != "2":
376	Divide(LHSWorkspace=nlist[i]+"pol",RHSWorkspace=P0nlist[i]+"pol",OutputWorkspace=nlist[i]+"SESANS")
377	if nspec > 4 and minSp != 3:
378	Divide(LHSWorkspace=nlist[i]+"2dpol",RHSWorkspace=P0nlist[i]+"pol",OutputWorkspace=nlist[i]+"2dSESANS")
379	else:
380	Divide(LHSWorkspace=nlist[i]+"pol",RHSWorkspace=P0nlist[i],OutputWorkspace=nlist[i]+"SESANS")
381	if nspec > 4 and minSp != 3:
382	Divide(LHSWorkspace=nlist[i]+"2dpol",RHSWorkspace=P0nlist[i],OutputWorkspace=nlist[i]+"2dSESANS")
383	ReplaceSpecialValues(InputWorkspace=nlist[i]+"SESANS",OutputWorkspace=nlist[i]+"SESANS",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
384
385	SEConstList=parseNameList(SEConstants)
386	k=0
387	for i in nlist:
388	if lnPOverLam == "2":
389	CloneWorkspace(InputWorkspace=i+"SESANS",OutputWorkspace=i+"SESANS_P")
390	a1=mtd[i+"SESANS"]
391	x=numpy.array(a1.readX(0))
392	new_y = numpy.array(a1.dataY(0))
393	new_e = numpy.array(a1.dataE(0))
394
395	for j in range(len(x)-1):
396	lam=((a1.readX(0)[j]+a1.readX(0)[j+1])/2.0)/10.0
397	p=a1.readY(0)[j]
398	e=a1.readE(0)[j]
399	if lnPOverLam == "2":
400	if p > 0.0:
401	new_y[j]=log(p)/((lam)**2)
402	new_e(0)[j]=(e/p)/((lam)**2)
403	else:
404	new_y(0)[j]=0.0
405	new_e(0)[j]=0.0
406	for j in range(len(x)):
407	if convertToSEL == "2":
408	lam=a1.readX(0)[j]
409	x[j]=1.0e-2float(SEConstList[k])lam*lam
410	#print str(lam)+" "+str(1.0e-2float(SEConstList[k])lam*lam)
411	a1.setY(0, new_y)
412	a1.setE(0, new_e)
413	a1.setX(0, x)
414	k=k+1
415
416
417	if nspec > 4 and minSp != 3:
418	k=0
419	for i in nlist:
420	if lnPOverLam == "2":
421	CloneWorkspace(InputWorkspace=i+"2dSESANS",OutputWorkspace=i+"2dSESANS_P")
422	a1=mtd[i+"2dSESANS"]
423	nspec=a1.getNumberHistograms()
424	for l in range(nspec):
425	x = numpy.array(a1.readX(l))
426	new_y = numpy.array(a1.readY(l))
427	new_e = numpy.array(a1.readE(l))
428	for j in range(len(x)-1):
429	lam=((a1.readX(l)[j]+a1.readX(l)[j+1])/2.0)/10.0
430	p=a1.readY(l)[j]
431	e=a1.readE(l)[j]
432	if lnPOverLam == "2":
433	if p > 0.0:
434	new_y[j]=log(p)/((lam1.0e-9)*2)
435	new_e[j]=(e/p)/((lam1.0e-9)*2)
436	else:
437	new_y[j]=0.0
438	new_e[j]=0.0
439	for j in range(len(x)):
440	if convertToSEL == "2":
441	lam=a1.readX(l)[j]
442	x[j]=1.0e-2float(SEConstList[k])lam*lam
443	#print str(lam)+" "+str(1.0e-2float(SEConstList[k])lam*lam)
444	a1.setY(l, new_y)
445	a1.setE(l, new_e)
446	a1.setX(l, x)
447	k=k+1
448
449	#
450	#===========================================================
451	#
452	def nrCalcSEConst(RFFrequency,poleShoeAngle):
453	if (RFFrequency=="0.5"):
454	B=0.53*34.288
455	elif (RFFrequency=="1.0"):
456	B=34.288
457	else:
458	B=2.0*34.288
459
460	h=6.62607e-34
461	m=1.67493e-27
462	L=1.0
463	Gl=1.83247e8
464	#
465	# correct the angle
466	# calibration of th0 using gold grating Dec 2010
467	#
468	th0=float(poleShoeAngle)
469	th0=-0.0000000467796(th05)+0.0000195413(th0*4)-0.00326229(th0*3)+0.271767(th0*2)-10.4269th0+198.108
470	c1=Glm2.0BL/(2.0pihtan(th0pi/180.0)*1.0e20)
471	print c1*1e8
472	return c1*1e8
473	#
474	#===========================================================
475	#
476	def nrSESANSP0Fn(P0runList,P0nameList,minSpec,maxSpec,upPeriod,downPeriod,gparams,diagnostics="0"):
477
478	P0nlist=parseNameList(P0nameList)
479	logger.notice("This is the P0nameslist:"+str(P0nlist))
480	P0rlist=parseRunList(P0runList)
481	logger.notice("This is the P0runlist:"+str(P0rlist))
482	for i in range(len(P0rlist)):
483	addRuns(P0rlist[i],P0nlist[i])
484
485	mon_spec=int(gparams[3])-1
486	minSp=int(minSpec)-1
487	maxSp=int(maxSpec)-1
488	reb=gparams[0]+","+gparams[1]+","+gparams[2]
489	if len(gparams) == 5:
490	mapfile=gparams[4]
491
492	for i in P0nlist:
493	a1=mtd[i+"_1"]
494	nspec=a1.getNumberHistograms()
495	if nspec == 1030 and minSp != 3:
496	GroupDetectors(InputWorkspace=i,OutputWorkspace=i,MapFile=mapfile)
497	ConvertUnits(InputWorkspace=i,OutputWorkspace=i,Target="Wavelength",AlignBins=1)
498	Rebin(InputWorkspace=i,OutputWorkspace=i,Params=reb)
499	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
500	if int(maxSpec) > int(minSpec):
501	SumSpectra(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=minSp,EndWorkspaceIndex=maxSp)
502	else:
503	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=minSp,EndWorkspaceIndex=maxSp)
504	Divide(LHSWorkspace=i+"det",RHSWorkspace=i+"mon",OutputWorkspace=i+"norm")
505	if (diagnostics=="0"):
506	DeleteWorkspace(i+"mon")
507	DeleteWorkspace(i+"det")
508	DeleteWorkspace(i)
509	Minus(LHSWorkspace=i+"norm_"+upPeriod,RHSWorkspace=i+"norm_"+downPeriod,OutputWorkspace="num")
510	Plus(LHSWorkspace=i+"norm_2",RHSWorkspace=i+"norm_1",OutputWorkspace="den")
511	Divide(LHSWorkspace="num",RHSWorkspace="den",OutputWorkspace=i+"pol")
512	ReplaceSpecialValues(InputWorkspace=i+"pol",OutputWorkspace=i+"pol",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
513	if (diagnostics=="0"):
514	DeleteWorkspace(i+"norm_2")
515	DeleteWorkspace(i+"norm_1")
516	DeleteWorkspace("num")
517	DeleteWorkspace("den")
518	#
519	#===========================================================
520	#
521	def nrSERGISFn(runList,nameList,P0runList,P0nameList,incidentAngles,SEConstants,specChan,minSpec,maxSpec,upPeriod,downPeriod,existingP0,gparams,lnPOverLam):
522	nlist=parseNameList(nameList)
523	logger.notice("This is the sample nameslist:"+str(nlist))
524	rlist=parseRunList(runList)
525	logger.notice("This is the sample runlist:"+str(rlist))
526	incAngles=parseNameList(incidentAngles)
527	logger.notice("This incident Angles are:"+str(incAngles))
528
529	for i in range(len(rlist)):
530	addRuns(rlist[i],nlist[i])
531
532	P0nlist=parseNameList(P0nameList)
533	logger.notice("This is the P0nameslist:"+str(P0nlist))
534	if existingP0 != "2":
535	P0rlist=parseRunList(P0runList)
536	logger.notice("This is the P0runlist:"+str(P0rlist))
537	for i in range(len(P0rlist)):
538	addRuns(P0rlist[i],P0nlist[i])
539
540	mon_spec=int(gparams[3])-1
541	minSp=int(minSpec)-1
542	maxSp=int(maxSpec)-1
543	reb=gparams[0]+","+gparams[1]+","+gparams[2]
544
545	k=0
546	for i in nlist:
547	for j in range(2):
548	wksp=i+"_"+pnums[j]
549	ConvertUnits(InputWorkspace=wksp,OutputWorkspace=wksp,Target="Wavelength",AlignBins=1)
550	Rebin(InputWorkspace=wksp,OutputWorkspace=wksp,Params=reb)
551	CropWorkspace(InputWorkspace=wksp,OutputWorkspace=wksp+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
552	a1=mtd[wksp]
553	nspec=a1.getNumberHistograms()
554	if nspec == 4:
555	CropWorkspace(InputWorkspace=wksp,OutputWorkspace=wksp+"det",StartWorkspaceIndex=3,EndWorkspaceIndex=3)
556	RotateInstrumentComponent(wksp+"det","DetectorBench",X="-1.0",Angle=str(2.0*float(incAngles[k])))
557	Divide(LHSWorkspace=wksp+"det",RHSWorkspace=wksp+"mon",OutputWorkspace=wksp+"norm")
558	else:
559	CropWorkspace(InputWorkspace=wksp,OutputWorkspace=wksp+"det",StartWorkspaceIndex=4,EndWorkspaceIndex=243)
560	# move the first spectrum in the list onto the beam centre so that when the bench is rotated it's in the right place
561	MoveInstrumentComponent(wksp+"det","DetectorBench",Y=str((125.0-float(minSpec))*1.2e-3))
562	# add a bit to the angle to put the first spectrum of the group in the right place
563	a1=2.0float(incAngles[k])+atan((float(minSpec)-float(specChan))1.2e-3/3.53)*180.0/pi
564	#print str(2.0float(incAngles[k]))+" "+str(atan((float(minSpec)-float(specChan))1.2e-3/3.63)*180.0/pi)+" "+str(a1)
565	RotateInstrumentComponent(wksp+"det","DetectorBench",X="-1.0",Angle=str(a1))
566	GroupDetectors(InputWorkspace=wksp+"det",OutputWorkspace=wksp+"sum",WorkspaceIndexList=range(int(minSpec)-5,int(maxSpec)-5+1),KeepUngroupedSpectra="0")
567	Divide(LHSWorkspace=wksp+"sum",RHSWorkspace=wksp+"mon",OutputWorkspace=wksp+"norm")
568	Divide(LHSWorkspace=wksp+"det",RHSWorkspace=wksp+"mon",OutputWorkspace=wksp+"detnorm")
569	floodnorm(wksp+"detnorm",floodfile)
570	DeleteWorkspace(wksp+"sum")
571
572	DeleteWorkspace(wksp+"mon")
573	DeleteWorkspace(wksp+"det")
574	DeleteWorkspace(wksp)
575
576	Minus(LHSWorkspace=i+"_"+upPeriod+"norm",RHSWorkspace=i+"_"+downPeriod+"norm",OutputWorkspace="num")
577	Plus(LHSWorkspace=i+"_1norm",RHSWorkspace=i+"_2norm",OutputWorkspace="den")
578	Divide(LHSWorkspace="num",RHSWorkspace="den",OutputWorkspace=i+"pol")
579	ReplaceSpecialValues(InputWorkspace=i+"pol",OutputWorkspace=i+"pol",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
580	DeleteWorkspace(i+"_1norm")
581	DeleteWorkspace(i+"_2norm")
582	DeleteWorkspace("num")
583	DeleteWorkspace("den")
584
585	if nspec != 4:
586	Minus(LHSWorkspace=i+"_"+upPeriod+"detnorm",RHSWorkspace=i+"_"+downPeriod+"detnorm",OutputWorkspace="num")
587	Plus(LHSWorkspace=i+"_1detnorm",RHSWorkspace=i+"_2detnorm",OutputWorkspace="den")
588	Divide(LHSWorkspace="num",RHSWorkspace="den",OutputWorkspace=i+"2dpol")
589	ReplaceSpecialValues(InputWorkspace=i+"2dpol",OutputWorkspace=i+"2dpol",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
590	DeleteWorkspace(i+"_1detnorm")
591	DeleteWorkspace(i+"_2detnorm")
592	DeleteWorkspace("num")
593	DeleteWorkspace("den")
594
595	if existingP0 != "2":
596	for i in P0nlist:
597	ConvertUnits(InputWorkspace=i,OutputWorkspace=i,Target="Wavelength",AlignBins=1)
598	Rebin(InputWorkspace=i,OutputWorkspace=i,Params=reb)
599	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
600	if int(maxSpec) > int(minSpec):
601	SumSpectra(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=minSp,EndWorkspaceIndex=maxSp)
602	else:
603	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=minSp,EndWorkspaceIndex=maxSp)
604	Divide(LHSWorkspace=i+"det",RHSWorkspace=i+"mon",OutputWorkspace=i+"norm")
605	DeleteWorkspace(i+"mon")
606	DeleteWorkspace(i+"det")
607	DeleteWorkspace(i)
608	Minus(LHSWorkspace=i+"norm_"+upPeriod,RHSWorkspace=i+"norm_"+downPeriod,OutputWorkspace="num")
609	Plus(LHSWorkspace=i+"norm_2",RHSWorkspace=i+"norm_1",OutputWorkspace="den")
610	Divide(LHSWorkspace="num",RHSWorkspace="den",OutputWorkspace=i+"pol")
611	ReplaceSpecialValues(InputWorkspace=i+"pol",OutputWorkspace=i+"pol",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
612	DeleteWorkspace(i+"norm_2")
613	DeleteWorkspace(i+"norm_1")
614	DeleteWorkspace("num")
615	DeleteWorkspace("den")
616
617	for i in range(len(nlist)):
618	if existingP0 != "2":
619	Divide(LHSWorkspace=nlist[i]+"pol",RHSWorkspace=P0nlist[i]+"pol",OutputWorkspace=nlist[i]+"SESANS")
620	else:
621	Divide(LHSWorkspace=nlist[i]+"pol",RHSWorkspace=P0nlist[i]+"pol",OutputWorkspace=nlist[i]+"SESANS")
622	ReplaceSpecialValues(InputWorkspace=nlist[i]+"SESANS",OutputWorkspace=nlist[i]+"SESANS",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
623
624	SEConstList=parseNameList(SEConstants)
625	k=0
626	for i in nlist:
627	a1=mtd[i+"SESANS"]
628	lam = ((a1.readX(0)[1:] + a1.readX(0)[:-1])/2.0)/10.0
629	p = a1.readY(0)
630	a1.setY(0, numpy.log(p)/((lam1.0e-8)*2) )
631	a1.setX(0, 1.0e-2float(SEConstList[k])lam*lam)
632	k=k+1
633	#
634	#===========================================================
635	#
636	def nrNRFn(runList,nameList,incidentAngles,DBList,specChan,minSpec,maxSpec,gparams,floodfile,subbgd=0,qgroup=0,dofloodnorm=1,diagnostics=0):
637	nlist=parseNameList(nameList)
638	#logger.notice("This is the sample nameslist:"+str(nlist))
639	rlist=parseRunList(runList)
640	#logger.notice("This is the sample runlist:"+str(rlist))
641	dlist=parseNameList(DBList)
642	#logger.notice("This is the Direct Beam nameslist:"+str(dlist))
643	incAngles=parseNameList(incidentAngles)
644	#logger.notice("This incident Angles are:"+str(incAngles))
645
646	for i in range(len(rlist)):
647	addRuns(rlist[i],nlist[i])
648
649	mon_spec=int(gparams[3])-1
650	reb=gparams[0]+","+gparams[1]+","+gparams[2]
651
652	k=0
653	for i in nlist:
654	if isinstance(mtd[i], WorkspaceGroup):
655	#RenameWorkspace(i+"_1",i)
656	snames=mtd[i].getNames()
657	Plus(LHSWorkspace=i+"_1",RHSWorkspace=i+"_2",OutputWorkspace="wtemp")
658	if len(snames) > 2:
659	for j in range(2,len(snames)-1):
660	Plus(LHSWorkspace="wtemp",RHSWorkspace=snames[j],OutputWorkspace="wtemp")
661	for j in snames:
662	DeleteWorkspace(j)
663	RenameWorkspace(InputWorkspace="wtemp",OutputWorkspace=i)
664	ConvertUnits(InputWorkspace=i,OutputWorkspace=i,Target="Wavelength",AlignBins="1")
665	a1=mtd[i]
666	nspec=a1.getNumberHistograms()
667	if nspec == 4:
668	Rebin(InputWorkspace=i,OutputWorkspace=i,Params=reb)
669	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
670	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=3,EndWorkspaceIndex=3)
671	RotateInstrumentComponent(i+"det","DetectorBench",X="-1.0",Angle=str(2.0*float(incAngles[k])))
672	Divide(LHSWorkspace=i+"det",RHSWorkspace=i+"mon",OutputWorkspace=i+"norm")
673	if dlist[k] != "none":
674	Divide(LHSWorkspace=i+"norm",RHSWorkspace=dlist[k],OutputWorkspace=i+"norm")
675	ReplaceSpecialValues(InputWorkspace=i+"norm",OutputWorkspace=i+"norm",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
676	ConvertUnits(InputWorkspace=i+"norm",OutputWorkspace=i+"RvQ",Target="MomentumTransfer")
677	else:
678	# Rebin using internal parameters to avoid problems with summing in Q
679	#internalreb=gparams[0]+",0.01,"+gparams[2]
680	#Rebin(InputWorkspace=i,OutputWorkspace=i,Params=internalreb)
681	minSp=int(minSpec)
682	maxSp=int(maxSpec)
683	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=4,EndWorkspaceIndex=243)
684	if(dofloodnorm==1):
685	floodnorm(i+"det",floodfile)
686	# move the first spectrum in the list onto the beam centre so that when the bench is rotated it's in the right place
687	MoveInstrumentComponent(i+"det","DetectorBench",Y=str((125.0-float(minSpec))*1.2e-3))
688	# add a bit to the angle to put the first spectrum of the group in the right place
689	a1=2.0float(incAngles[k])+atan((float(minSpec)-float(specChan))1.2e-3/3.53)*180.0/pi
690	#print str(2.0float(incAngles[k]))+" "+str(atan((float(minSpec)-float(specChan))1.2e-3/3.63)*180.0/pi)+" "+str(a1)
691	RotateInstrumentComponent(i+"det","DetectorBench",X="-1.0",Angle=str(a1))
692	if (subbgd==1):
693	# Calculate a background correction
694	GroupDetectors(i+"det",OutputWorkspace=i+"bgd2",WorkspaceIndexList=range(0,50),KeepUngroupedSpectra="0")
695	GroupDetectors(i+"det",OutputWorkspace=i+"bgd1",WorkspaceIndexList=range(160,240),KeepUngroupedSpectra="0")
696	Plus(i+"bgd1",i+"bgd2",OutputWorkspace=i+"bgd")
697	wbgdtemp=mtd[i+"bgd"]/130.0
698	DeleteWorkspace(i+"bgd1")
699	DeleteWorkspace(i+"bgd2")
700	DeleteWorkspace(i+"bgd")
701	# Subract a per spectrum background
702	Minus(i+"det",wbgdtemp,OutputWorkspace=i+"det")
703	if(diagnostics==0):
704	DeleteWorkspace("wbgdtemp")
705	# Experimental convert to Q before summing
706	if (qgroup==1):
707	Rebin(InputWorkspace=i+"det",OutputWorkspace=i+"det",Params=reb)
708	CropWorkspace(InputWorkspace=i+"det",OutputWorkspace=i+"detQ",StartWorkspaceIndex=int(minSpec)-5,EndWorkspaceIndex=int(maxSpec)-5+1)
709	ConvertUnits(InputWorkspace=i+"detQ",OutputWorkspace=i+"detQ",Target="MomentumTransfer",AlignBins="1")
710	GroupDetectors(i+"detQ",OutputWorkspace=i+"sum",WorkspaceIndexList=range(0,int(maxSpec)-int(minSpec)),KeepUngroupedSpectra="0")
711	ConvertUnits(InputWorkspace=i+"sum",OutputWorkspace=i+"sum",Target="Wavelength",AlignBins="1")
712	Rebin(InputWorkspace=i+"sum",OutputWorkspace=i+"sum",Params=reb)
713	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
714	Rebin(InputWorkspace=i+"mon",OutputWorkspace=i+"mon",Params=reb)
715	Rebin(InputWorkspace=i+"det",OutputWorkspace=i+"det",Params=reb)
716	DeleteWorkspace(i+"detQ")
717	else:
718	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
719	Rebin(InputWorkspace=i+"mon",OutputWorkspace=i+"mon",Params=reb)
720	Rebin(InputWorkspace=i+"det",OutputWorkspace=i+"det",Params=reb)
721	GroupDetectors(InputWorkspace=i+"det",OutputWorkspace=i+"sum",WorkspaceIndexList=range(int(minSpec)-5,int(maxSpec)-5+1),KeepUngroupedSpectra="0")
722	Divide(LHSWorkspace=i+"sum",RHSWorkspace=i+"mon",OutputWorkspace=i+"norm")
723	Divide(LHSWorkspace=i+"det",RHSWorkspace=i+"mon",OutputWorkspace=i+"detnorm")
724	if dlist[k] == "none":
725	a1=0
726	elif dlist[k] == "function":
727	# polynomial + power law corrections based on Run numbers 8291 and 8292
728	Divide(LHSWorkspace=i+'norm',RHSWorkspace=i+'norm',OutputWorkspace=i+'normt1')
729	PolynomialCorrection(InputWorkspace=i+'normt1',OutputWorkspace=i+'normPC',Coefficients='-0.0177398,0.00101695,0.0',Operation='Multiply')
730	PowerLawCorrection(InputWorkspace=i+'normt1',OutputWorkspace=i+'normPLC',C0='2.01332',C1='-1.8188')
731	Plus(LHSWorkspace=i+'normPC',RHSWorkspace=i+'normPLC',OutputWorkspace=i+'normt1')
732	Divide(LHSWorkspace=i+'norm',RHSWorkspace=i+'normt1',OutputWorkspace=i+'norm')
733	ReplaceSpecialValues(InputWorkspace=i+'norm',OutputWorkspace=i+'norm',NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
734	DeleteWorkspace(i+'normPC')
735	DeleteWorkspace(i+'normPLC')
736	DeleteWorkspace(i+'normt1')
737	else:
738	Divide(LHSWorkspace=i+"norm",RHSWorkspace=dlist[k],OutputWorkspace=i+"norm")
739	ReplaceSpecialValues(InputWorkspace=i+"norm",OutputWorkspace=i+"norm",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
740	Divide(LHSWorkspace=i+"detnorm",RHSWorkspace=dlist[k],OutputWorkspace=i+"detnorm")
741	ReplaceSpecialValues(InputWorkspace=i+"detnorm",OutputWorkspace=i+"detnorm",NanValue=0.0,NaNError=0.0,InfinityValue=0.0,InfinityError=0.0)
742	ConvertUnits(InputWorkspace=i+"norm",OutputWorkspace=i+"RvQ",Target="MomentumTransfer")
743	#floodnorm(i+"detnorm",floodfile)
744	DeleteWorkspace(i+"sum")
745
746	k=k+1
747	DeleteWorkspace(i)
748	if(diagnostics==0):
749	DeleteWorkspace(i+"mon")
750	DeleteWorkspace(i+"det")
751
752	#
753	#===========================================================
754	#
755	def findbin(wksp,val):
756	a1=mtd[wksp]
757	x1=a1.readX(0)
758	bnum=-1
759	for i in range(len(x1)-1):
760	if x1[i] > val:
761	break
762	return i-1
763	#
764	#===========================================================
765	#
766	def nrDBFn(runListShort,nameListShort,runListLong,nameListLong,nameListComb,minSpec,maxSpec,minWavelength,gparams,floodfile="",diagnostics="0"):
767	nlistS=parseNameList(nameListShort)
768	rlistS=parseRunList(runListShort)
769	nlistL=parseNameList(nameListLong)
770	rlistL=parseRunList(runListLong)
771	nlistComb=parseNameList(nameListComb)
772
773	for i in range(len(rlistS)):
774	addRuns(rlistS[i],nlistS[i])
775	for i in range(len(rlistL)):
776	addRuns(rlistL[i],nlistL[i])
777
778	mon_spec=int(gparams[3])-1
779	minSp=int(minSpec)-1
780	maxSp=int(maxSpec)-1
781	reb=gparams[0]+","+gparams[1]+","+gparams[2]
782
783	for i in nlistS:
784	ConvertUnits(InputWorkspace=i,OutputWorkspace=i,Target="Wavelength",AlignBins="1")
785	Rebin(InputWorkspace=i,OutputWorkspace=i,Params=reb)
786	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
787	if isinstance(mtd[i], GroupWorkspace):
788	snames=mtd[i].getNames()
789	a1=mtd[snames[0]]
790	else:
791	a1=mtd[i]
792
793	nspec=a1.getNumberHistograms()
794
795	if nspec == 4:
796	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=3,EndWorkspaceIndex=3)
797	Divide(i+"det",i+"mon",OutputWorkspace=i+"norm")
798	ReplaceSpecialValues(i+"norm","0.0","0.0","0.0","0.0",OutputWorkspace=i+"norm")
799	else:
800	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=4,EndWorkspaceIndex=243)
801	floodnorm(i+"det",floodfile)
802	GroupDetectors(i+"det",OutputWorkspace=i+"sum",WorkspaceIndexList=range(int(minSpec)-5,int(maxSpec)-5+1),KeepUngroupedSpectra="0")
803	Divide(i+"sum",i+"mon",OutputWorkspace=i+"norm")
804	ReplaceSpecialValues(i+"norm","0.0","0.0","0.0","0.0", OutputWorkspace=i+"norm")
805
806	for i in nlistL:
807	ConvertUnits(InputWorkspace=i,OutputWorkspace=i,Target="Wavelength",AlignBins="1")
808	Rebin(InputWorkspace=i,OutputWorkspace=i,Params=reb)
809	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
810	if isinstance(mtd[i], GroupWorkspace):
811	lnames=mtd[i].getNames()
812	a1=mtd[lnames[0]]
813	else:
814	a1=mtd[i]
815
816	nspec=a1.getNumberHistograms()
817
818	if nspec == 4:
819	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=3,EndWorkspaceIndex=3)
820	Divide(i+"det",i+"mon",OutputWorkspace=i+"norm")
821	ReplaceSpecialValues(i+"norm","0.0","0.0","0.0","0.0",OutputWorkspace=i+"norm")
822	else:
823	CropWorkspace(InputWorkspace=i,OutputWorkspace=i+"det",StartWorkspaceIndex=4,EndWorkspaceIndex=243)
824	floodnorm(i+"det",floodfile)
825	GroupDetectors(i+"det",OutputWorkspace=i+"sum",WorkspaceIndexList=range(int(minSpec)-5,int(maxSpec)-5+1),KeepUngroupedSpectra="0")
826	Divide(i+"sum",i+"mon",OutputWorkspace=i+"norm")
827	ReplaceSpecialValues(i+"norm","0.0","0.0","0.0","0.0",OutputWorkspace=i+"norm")
828
829	for i in range(len(nlistS)):
830	if isinstance(mtd[nlistS[i]+"norm"], GroupWorkspace):
831	snames=mtd[nlistS[i]+"norm"].getNames()
832	lnames=mtd[nlistL[i]+"norm"].getNames()
833	for k in range(len(snames)):
834	Integration(snames[k],minWavelength,gparams[2], OutputWorkspace=snames[k]+"int")
835	Integration(lnames[k],minWavelength,gparams[2], OutputWorkspace=lnames[k]+"int")
836	Multiply(snames[k],lnames[k]+"int",OutputWorkspace=snames[k])
837	Divide(snames[k],snames[k]+"int",OutputWorkspace=snames[k])
838	a1=findbin(lnames[k],float(minWavelength))
839	MultiplyRange(lnames[k],"0",str(a1),"0.0",OutputWorkspace=lnames[k])
840	WeightedMean(snames[k],nlistComb[i]+"_"+str(k+1),OutputWorkspace=lnames[k])
841	if (diagnostics=="0"):
842	DeleteWorkspace(snames[k]+"int")
843	DeleteWorkspace(lnames[k]+"int")
844	else:
845	Integration(nlistS[i]+"norm",minWavelength,gparams[2],OutputWorkspace=nlistS[i]+"int")
846	Integration(nlistL[i]+"norm",minWavelength,gparams[2],OutputWorkspace=nlistL[i]+"int")
847	Multiply(nlistS[i]+"norm",nlistL[i]+"int",OutputWorkspace=nlistS[i]+"norm")
848	Divide(nlistS[i]+"norm",nlistS[i]+"int",OutputWorkspace=nlistS[i]+"norm")
849	a1=findbin(nlistL[i]+"norm",float(minWavelength))
850	MultiplyRange(nlistL[i]+"norm","0",str(a1),"0.0", OutputWorkspace=nlistL[i]+"norm")
851	WeightedMean(nlistS[i]+"norm",nlistL[i]+"norm", OutputWorkspace=nlistComb[i])
852	if (diagnostics=="0"):
853	DeleteWorkspace(nlistS[i]+"int")
854	DeleteWorkspace(nlistL[i]+"int")
855
856	if (diagnostics=="0"):
857	DeleteWorkspace(nlistS[i]+"mon")
858	DeleteWorkspace(nlistS[i]+"det")
859	if nspec != 4:
860	DeleteWorkspace(nlistS[i]+"sum")
861	DeleteWorkspace(nlistS[i]+"norm")
862	DeleteWorkspace(nlistS[i])
863	DeleteWorkspace(nlistL[i]+"mon")
864	DeleteWorkspace(nlistL[i]+"det")
865	if nspec != 4:
866	DeleteWorkspace(nlistL[i]+"sum")
867	DeleteWorkspace(nlistL[i]+"norm")
868	DeleteWorkspace(nlistL[i])
869
870	#
871	#===========================================================
872	#
873	def numberofbins(wksp):
874	a1=mtd[wksp]
875	y1=a1.readY(0)
876	return len(y1)-1
877	#
878	#===========================================================
879	#
880	def maskbin(wksp,val):
881	a1=mtd[wksp]
882	x1=a1.readX(0)
883	for i in range(len(x1)-1):
884	if x1[i] > val:
885	break
886	a1.dataY(0)[i-1]=0.0
887	a1.dataE(0)[i-1]=0.0
888	#
889	#===========================================================
890	#
891	def arr2list(iarray):
892	# convert array of strings to a single string with commas
893	res=""
894	for i in range(len(iarray)-1):
895	res=res+iarray[i]+","
896	res=res+iarray[len(iarray)-1]
897	return res
898	#
899	#===========================================================
900	#
901	def NRCombineDatafn(RunsNameList,CombNameList,applySFs,SFList,SFError,scaleOption,bparams,globalSF,applyGlobalSF,diagnostics=0):
902	qmin=bparams[0]
903	bin=bparams[1]
904	qmax=bparams[2]
905	rlist=parseNameList(RunsNameList)
906	listsfs=parseNameList(SFList)
907	listsfserr=parseNameList(SFError)
908	sfs=[]
909	sferrs=[]
910	for i in rlist:
911	Rebin(i,qmin+","+bin+","+qmax,OutputWorkspace=i+"reb")
912	# find the overlap ranges
913	bol=[] #beginning of overlaps
914	eol=[] #end of overlaps
915	for i in range(len(rlist)-1):
916	a1=mtd[rlist[i+1]]
917	x=a1.readX(0)
918	bol.append(x[0])
919	a1=mtd[rlist[i]]
920	x=a1.readX(0)
921	eol.append(x[len(x)-1])
922	# set the edges of the rebinned data to 0.0 to avoid partial bin problems
923	maskbin(rlist[0]+"reb",eol[0])
924	if len(rlist) > 2:
925	for i in range(1,len(rlist)-1):
926	maskbin(rlist[i]+"reb",bol[i-1])
927	maskbin(rlist[i]+"reb",eol[i])
928	maskbin(rlist[len(rlist)-1]+"reb",bol[len(rlist)-2])
929	# Now find the various scale factors and store in temp workspaces
930	for i in range(len(rlist)-1):
931	Integration(rlist[i]+"reb",str(bol[i]),str(eol[i]), OutputWorkspace="i"+str(i)+"1temp")
932	Integration(rlist[i+1]+"reb",str(bol[i]),str(eol[i]), OutputWorkspace="i"+str(i)+"2temp")
933	if scaleOption != "2":
934	Divide("i"+str(i)+"1temp","i"+str(i)+"2temp",OutputWorkspace="sf"+str(i))
935	a1=mtd["sf"+str(i)]
936	print "sf"+str(i)+"="+str(a1.readY(0))+" +/- "+str(a1.readE(0))
937	sfs.append(str(a1.readY(0)[0]))
938	sferrs.append(str(a1.readE(0)[0]))
939	else:
940	Divide("i"+str(i)+"2temp","i"+str(i)+"1temp",OutputWorkspace="sf"+str(i))
941	print "sf"+str(i)+"="+str(a1.readY(0))+" +/- "+str(a1.readE(0))
942	sfs.append(str(a1.readY(0)[0]))
943	sferrs.append(str(a1.readE(0)[0]))
944	DeleteWorkspace("i"+str(i)+"1temp")
945	DeleteWorkspace("i"+str(i)+"2temp")
946	# if applying pre-defined scale factors substitute the given values now
947	# Note the errors are now set to 0
948	if applySFs == "2":
949	for i in range(len(rlist)-1):
950	a1=mtd["sf"+str(i)]
951	a1.setY(0, float(listsfs[i]))
952	a1.setE(0, float(listsfserr[i]))
953	# Now scale the various data sets in the correct order
954	if scaleOption != "2":
955	for i in range(len(rlist)-1):
956	for j in range(i+1,len(rlist)):
957	Multiply(rlist[j]+"reb","sf"+str(i),OutputWorkspace=rlist[j]+"reb")
958	else:
959	for i in range(len(rlist)-1,0,-1):
960	for j in range(i,0,-1):
961	Multiply(rlist[j]+"reb","sf"+str(i-1),OutputWorkspace=rlist[j]+"reb")
962
963	WeightedMean(rlist[0]+"reb",rlist[1]+"reb",OutputWorkspace="currentSum")
964	if len(rlist) > 2:
965	for i in range(2,len(rlist)):
966	WeightedMean("currentSum",rlist[i]+"reb",OutputWorkspace="currentSum")
967
968	# if applying a global scale factor do it here
969	if applyGlobalSF == "2":
970	scaledData=mtd['currentSum']/float(globalSF)
971	RenameWorkspace('scaledData',OutputWorkspace=CombNameList)
972	DeleteWorkspace('currentSum')
973	else:
974	RenameWorkspace('currentSum',OutputWorkspace=CombNameList)
975	for i in range(len(rlist)-1):
976	DeleteWorkspace("sf"+str(i))
977	if (diagnostics==0):
978	for i in range(len(rlist)):
979	DeleteWorkspace(rlist[i]+"reb")
980	return [arr2list(sfs),arr2list(sferrs)]
981
982	#
983	#===========================================================
984	#
985	def nrWriteXYE(wksp,fname):
986	a1=mtd[wksp]
987	x1=a1.readX(0)
988	X1=n.zeros((len(x1)-1))
989	for i in range(0,len(x1)-1):
990	X1[i]=(x1[i]+x1[i+1])/2.0
991	y1=a1.readY(0)
992	e1=a1.readE(0)
993	f=open(fname,'w')
994	for i in range(len(X1)):
995	s=""
996	s+="%f," % X1[i]
997	s+="%f," % y1[i]
998	s+="%f\n" % e1[i]
999	f.write(s)
1000	f.close()
1001	#
1002	#===========================================================
1003	#
1004	def nrPNRCorrection(UpWksp,DownWksp,calibration=0):
1005	# crho=[0.941893,0.0234006,-0.00210536,0.0]
1006	# calpha=[0.945088,0.0242861,-0.00213624,0.0]
1007	# cAp=[1.00079,-0.0186778,0.00131546,0.0]
1008	# cPp=[1.01649,-0.0228172,0.00214626,0.0]
1009	# Constants Based on Runs 18350+18355 and 18351+18356 analyser theta at -0.1deg
1010	# 2 RF Flippers as the polarising system
1011	if (calibration == 0):
1012	crho=[1.006831,-0.011467,0.002244,-0.000095]
1013	calpha=[1.017526,-0.017183,0.003136,-0.000140]
1014	cAp=[0.917940,0.038265,-0.006645,0.000282]
1015	cPp=[0.972762,0.001828,-0.000261,0.0]
1016	elif (calibration == 1):
1017	# Constants Based on Runs 19438-19458 and 19439+19459
1018	# Drabkin on incident side RF Flipper on analyser side, Dec 2012
1019	crho=[0.970257,0.016127,-0.002318,0.000090]
1020	calpha=[0.975722,0.012464,-0.002408,0.000105]
1021	cAp=[1.030894,-0.040847,0.006069,-0.000247]
1022	cPp=[0.961900,-0.003722,0.001094,-0.000057]
1023	elif (calibration == 2):
1024	# Constants Based on Runs 19628-19656 and 19660-19670 analyser -0.2deg
1025	# RF Flippers on polariser and analyser side Feb 2013
1026	crho=[0.945927,0.025421,-0.003647,0.000156]
1027	calpha=[0.940769,0.027250,-0.003848,0.000164]
1028	cAp=[0.974374,-0.005334,0.001313,-0.000115]
1029	cPp=[1.023141,-0.024548,0.003398,-0.000134]
1030	elif (calibration == 3):
1031	# Constants Based on Runs 19628-19656 and 19660-19670 analyser -0.1deg
1032	# RF Flippers on polariser and analyser side Feb 2013
1033	crho=[0.955384,0.021501,-0.002962,0.000112]
1034	calpha=[0.957789,0.019995,-0.002697,0.000099]
1035	cAp=[0.986906,-0.013945,0.002480,-0.000161]
1036	cPp=[0.999517,-0.013878,0.001680,-0.000043]
1037
1038	Ip = mtd[UpWksp]
1039	Ia = mtd[DownWksp]
1040	CloneWorkspace(Ip,OutputWorkspace="PCalpha")
1041	CropWorkspace(InputWorkspace="PCalpha",OutputWorkspace="PCalpha",StartWorkspaceIndex="0",EndWorkspaceIndex="0")
1042	PCalpha=(mtd['PCalpha']*0.0)+1.0
1043	alpha=mtd['PCalpha']
1044	# a1=alpha.readY(0)
1045	# for i in range(0,len(a1)):
1046	# alpha.dataY(0)[i]=0.0
1047	# alpha.dataE(0)[i]=0.0
1048	CloneWorkspace("PCalpha",OutputWorkspace="PCrho")
1049	CloneWorkspace("PCalpha",OutputWorkspace="PCAp")
1050	CloneWorkspace("PCalpha",OutputWorkspace="PCPp")
1051	rho=mtd['PCrho']
1052	Ap=mtd['PCAp']
1053	Pp=mtd['PCPp']
1054	# for i in range(0,len(a1)):
1055	# x=(alpha.dataX(0)[i]+alpha.dataX(0)[i])/2.0
1056	# for j in range(0,4):
1057	# alpha.dataY(0)[i]=alpha.dataY(0)[i]+calpha[j]x*j
1058	# rho.dataY(0)[i]=rho.dataY(0)[i]+crho[j]x*j
1059	# Ap.dataY(0)[i]=Ap.dataY(0)[i]+cAp[j]x*j
1060	# Pp.dataY(0)[i]=Pp.dataY(0)[i]+cPp[j]x*j
1061	PolynomialCorrection(InputWorkspace="PCalpha",OutputWorkspace="PCalpha",Coefficients=calpha,Operation="Multiply")
1062	PolynomialCorrection(InputWorkspace="PCrho",OutputWorkspace="PCrho",Coefficients=crho,Operation="Multiply")
1063	PolynomialCorrection(InputWorkspace="PCAp",OutputWorkspace="PCAp",Coefficients=cAp,Operation="Multiply")
1064	PolynomialCorrection(InputWorkspace="PCPp",OutputWorkspace="PCPp",Coefficients=cPp,Operation="Multiply")
1065	D=Pp*(1.0+rho)
1066	nIp=(Ip(rhoPp+1.0)+Ia*(Pp-1.0))/D
1067	nIa=(Ip(rhoPp-1.0)+Ia*(Pp+1.0))/D
1068	RenameWorkspace(nIp,OutputWorkspace=str(Ip)+"corr")
1069	RenameWorkspace(nIa,OutputWorkspace=str(Ia)+"corr")
1070	iwksp=mtd.getObjectNames()
1071	list_n=[str(Ip),str(Ia),"PCalpha","PCrho","PCAp","PCPp","1_p"]
1072	for i in range(len(iwksp)):
1073	for j in list_n:
1074	lname=len(j)
1075	if iwksp[i] [0:lname+1] == j+"_":
1076	DeleteWorkspace(iwksp[i])
1077	DeleteWorkspace("PCalpha")
1078	DeleteWorkspace("PCrho")
1079	DeleteWorkspace("PCAp")
1080	DeleteWorkspace("PCPp")
1081	DeleteWorkspace("D")
1082	#
1083	#===========================================================
1084	#
1085	def nrPACorrection(UpUpWksp,UpDownWksp,DownUpWksp,DownDownWksp,calibration=0):
1086	# crho=[0.941893,0.0234006,-0.00210536,0.0]
1087	# calpha=[0.945088,0.0242861,-0.00213624,0.0]
1088	# cAp=[1.00079,-0.0186778,0.00131546,0.0]
1089	# cPp=[1.01649,-0.0228172,0.00214626,0.0]
1090	# Constants Based on Runs 18350+18355 and 18351+18356 analyser theta at -0.1deg
1091	# 2 RF Flippers as the polarising system
1092	if (calibration == 0):
1093	crho=[1.006831,-0.011467,0.002244,-0.000095]
1094	calpha=[1.017526,-0.017183,0.003136,-0.000140]
1095	cAp=[0.917940,0.038265,-0.006645,0.000282]
1096	cPp=[0.972762,0.001828,-0.000261,0.0]
1097	elif (calibration == 1):
1098	# Constants Based on Runs 19438-19458 and 19439+19459
1099	# Drabkin on incident side RF Flipper on analyser side Dec 2012
1100	crho=[0.970257,0.016127,-0.002318,0.000090]
1101	calpha=[0.975722,0.012464,-0.002408,0.000105]
1102	cAp=[1.030894,-0.040847,0.006069,-0.000247]
1103	cPp=[0.961900,-0.003722,0.001094,-0.000057]
1104	elif (calibration == 2):
1105	# Constants Based on Runs 19628-19656 and 19660-19670 analyser -0.2deg
1106	# RF Flippers on polariser and analyser side Feb 2013
1107	crho=[0.945927,0.025421,-0.003647,0.000156]
1108	calpha=[0.940769,0.027250,-0.003848,0.000164]
1109	cAp=[0.974374,-0.005334,0.001313,-0.000115]
1110	cPp=[1.023141,-0.024548,0.003398,-0.000134]
1111	elif (calibration == 3):
1112	# Constants Based on Runs 19628-19656 and 19660-19670 analyser -0.1deg
1113	# RF Flippers on polariser and analyser side Feb 2013
1114	crho=[0.955384,0.021501,-0.002962,0.000112]
1115	calpha=[0.957789,0.019995,-0.002697,0.000099]
1116	cAp=[0.986906,-0.013945,0.002480,-0.000161]
1117	cPp=[0.999517,-0.013878,0.001680,-0.000043]
1118
1119	Ipp = mtd[UpUpWksp]
1120	Ipa = mtd[UpDownWksp]
1121	Iap = mtd[DownUpWksp]
1122	Iaa = mtd[DownDownWksp]
1123	CloneWorkspace(Ipp,OutputWorkspace="PCalpha")
1124	CropWorkspace(InputWorkspace="PCalpha",OutputWorkspace="PCalpha",StartWorkspaceIndex="0",EndWorkspaceIndex="0")
1125	PCalpha=(mtd['PCalpha']*0.0)+1.0
1126	alpha=mtd['PCalpha']
1127	# a1=alpha.readY(0)
1128	# for i in range(0,len(a1)):
1129	# alpha.dataY(0)[i]=0.0
1130	# alpha.dataE(0)[i]=0.0
1131	CloneWorkspace("PCalpha",OutputWorkspace="PCrho")
1132	CloneWorkspace("PCalpha",OutputWorkspace="PCAp")
1133	CloneWorkspace("PCalpha",OutputWorkspace="PCPp")
1134	rho=mtd['PCrho']
1135	Ap=mtd['PCAp']
1136	Pp=mtd['PCPp']
1137	# for i in range(0,len(a1)):
1138	# x=(alpha.dataX(0)[i]+alpha.dataX(0)[i])/2.0
1139	# for j in range(0,4):
1140	# alpha.dataY(0)[i]=alpha.dataY(0)[i]+calpha[j]x*j
1141	# rho.dataY(0)[i]=rho.dataY(0)[i]+crho[j]x*j
1142	# Ap.dataY(0)[i]=Ap.dataY(0)[i]+cAp[j]x*j
1143	# Pp.dataY(0)[i]=Pp.dataY(0)[i]+cPp[j]x*j
1144	# Use the polynomial corretion fn instead
1145	PolynomialCorrection(InputWorkspace="PCalpha",OutputWorkspace="PCalpha",Coefficients=calpha,Operation="Multiply")
1146	PolynomialCorrection(InputWorkspace="PCrho",OutputWorkspace="PCrho",Coefficients=crho,Operation="Multiply")
1147	PolynomialCorrection(InputWorkspace="PCAp",OutputWorkspace="PCAp",Coefficients=cAp,Operation="Multiply")
1148	PolynomialCorrection(InputWorkspace="PCPp",OutputWorkspace="PCPp",Coefficients=cPp,Operation="Multiply")
1149
1150	A0 = (Iaa * Pp * Ap) + (Ap * Ipa * rho * Pp) + (Ap * Iap * Pp * alpha) + (Ipp * Ap * alpha * rho * Pp)
1151	A1 = Pp * Iaa
1152	A2 = Pp * Iap
1153	A3 = Ap * Iaa
1154	A4 = Ap * Ipa
1155	A5 = Ap * alpha * Ipp
1156	A6 = Ap * alpha * Iap
1157	A7 = Pp * rho * Ipp
1158	A8 = Pp * rho * Ipa
1159	D = Pp * Ap ( 1.0 + rho + alpha + (rho alpha) )
1160	nIpp = (A0 - A1 + A2 - A3 + A4 + A5 - A6 + A7 - A8 + Ipp + Iaa - Ipa - Iap) / D
1161	nIaa = (A0 + A1 - A2 + A3 - A4 - A5 + A6 - A7 + A8 + Ipp + Iaa - Ipa - Iap) / D
1162	nIpa = (A0 - A1 + A2 + A3 - A4 - A5 + A6 + A7 - A8 - Ipp - Iaa + Ipa + Iap) / D
1163	nIap = (A0 + A1 - A2 - A3 + A4 + A5 - A6 - A7 + A8 - Ipp - Iaa + Ipa + Iap) / D
1164	RenameWorkspace(nIpp,OutputWorkspace=str(Ipp)+"corr")
1165	RenameWorkspace(nIpa,OutputWorkspace=str(Ipa)+"corr")
1166	RenameWorkspace(nIap,OutputWorkspace=str(Iap)+"corr")
1167	RenameWorkspace(nIaa,OutputWorkspace=str(Iaa)+"corr")
1168	ReplaceSpecialValues(str(Ipp)+"corr",OutputWorkspace=str(Ipp)+"corr",NaNValue="0.0",NaNError="0.0",InfinityValue="0.0",InfinityError="0.0")
1169	ReplaceSpecialValues(str(Ipp)+"corr",OutputWorkspace=str(Ipp)+"corr",NaNValue="0.0",NaNError="0.0",InfinityValue="0.0",InfinityError="0.0")
1170	ReplaceSpecialValues(str(Ipp)+"corr",OutputWorkspace=str(Ipp)+"corr",NaNValue="0.0",NaNError="0.0",InfinityValue="0.0",InfinityError="0.0")
1171	ReplaceSpecialValues(str(Ipp)+"corr",OutputWorkspace=str(Ipp)+"corr",NaNValue="0.0",NaNError="0.0",InfinityValue="0.0",InfinityError="0.0")
1172	iwksp=mtd.getObjectNames()
1173	list_n=[str(Ipp),str(Ipa),str(Iap),str(Iaa),"PCalpha","PCrho","PCAp","PCPp","1_p"]
1174	for i in range(len(iwksp)):
1175	for j in list_n:
1176	lname=len(j)
1177	if iwksp[i] [0:lname+1] == j+"_":
1178	DeleteWorkspace(iwksp[i])
1179	DeleteWorkspace("PCalpha")
1180	DeleteWorkspace("PCrho")
1181	DeleteWorkspace("PCAp")
1182	DeleteWorkspace("PCPp")
1183	DeleteWorkspace('A0')
1184	DeleteWorkspace('A1')
1185	DeleteWorkspace('A2')
1186	DeleteWorkspace('A3')
1187	DeleteWorkspace('A4')
1188	DeleteWorkspace('A5')
1189	DeleteWorkspace('A6')
1190	DeleteWorkspace('A7')
1191	DeleteWorkspace('A8')
1192	DeleteWorkspace('D')
1193	#
1194	#===========================================================
1195	#
1196	def nrPNRFn(runList,nameList,incidentAngles,DBList,specChan,minSpec,maxSpec,gparams,floodfile,PNRwithPA,pnums,doCorrs,doLDCorrs="0",subbgd=0,calibration=0,diagnostics=0):
1197	nlist=parseNameList(nameList)
1198	logger.notice("This is the sample nameslist:"+str(nlist))
1199	rlist=parseRunList(runList)
1200	logger.notice("This is the sample runlist:"+str(rlist))
1201	dlist=parseNameList(DBList)
1202	logger.notice("This is the Direct Beam nameslist:"+str(dlist))
1203	incAngles=parseNameList(incidentAngles)
1204	logger.notice("This incident Angles are:"+str(incAngles))
1205
1206	if PNRwithPA == "2":
1207	nper=4
1208	logger.notice("PNRwithPA="+str(PNRwithPA))
1209	logger.notice(str(pnums))
1210	else:
1211	nper=2
1212
1213	for i in range(len(rlist)):
1214	addRuns(rlist[i],nlist[i])
1215
1216	mon_spec=int(gparams[3])-1
1217	minSp=int(minSpec)
1218	maxSp=int(maxSpec)
1219	reb=gparams[0]+","+gparams[1]+","+gparams[2]
1220
1221	k=0
1222	for i in nlist:
1223	a1=mtd[i+"_1"]
1224	nspec=a1.getNumberHistograms()
1225	if (subbgd==1 and nspec!=4):
1226	# If a background subtraction is required sum the bgd outside the
1227	# area of the detector that is visible through the analyser over all periods and average
1228	CloneWorkspace(i, OutputWorkspace="bgdtemp")
1229	ConvertUnits(InputWorkspace="bgdtemp",OutputWorkspace="bgdtemp",Target="Wavelength",AlignBins="1")
1230	Rebin(InputWorkspace="bgdtemp",OutputWorkspace="bgdtemp",Params=reb)
1231	CropWorkspace(InputWorkspace="bgdtemp",OutputWorkspace="bgdtemp",StartWorkspaceIndex=4,EndWorkspaceIndex=243)
1232	Plus("bgdtemp"+"_"+pnums[0],"bgdtemp"+"_"+pnums[1],OutputWorkspace="wbgdsum")
1233	if (nper>2):
1234	for j in range(2,nper):
1235	Plus("wbgdsum","bgdtemp"+"_"+pnums[j],OutputWorkspace="wbgdsum")
1236	GroupDetectors("wbgdsum",OutputWorkspace="bgd2",WorkspaceIndexList=range(0,50),KeepUngroupedSpectra="0")
1237	GroupDetectors("wbgdsum",OutputWorkspace="bgd1",WorkspaceIndexList=range(160,240),KeepUngroupedSpectra="0")
1238	Plus("bgd1","bgd2",OutputWorkspace="bgd")
1239	wbgdtemp=mtd["bgd"]/(130.0*nper)
1240	DeleteWorkspace("bgdtemp")
1241	DeleteWorkspace("wbgdsum")
1242	DeleteWorkspace("bgd1")
1243	DeleteWorkspace("bgd2")
1244	DeleteWorkspace("bgd")
1245
1246	wksp=i
1247	ConvertUnits(InputWorkspace=wksp,OutputWorkspace=wksp,Target="Wavelength",AlignBins="1")
1248	Rebin(InputWorkspace=wksp,OutputWorkspace=wksp,Params=reb)
1249	CropWorkspace(InputWorkspace=wksp,OutputWorkspace=wksp+"mon",StartWorkspaceIndex=mon_spec,EndWorkspaceIndex=mon_spec)
1250	if nspec == 4:
1251	CropWorkspace(InputWorkspace=wksp,OutputWorkspace=wksp+"det",StartWorkspaceIndex=3,EndWorkspaceIndex=3)
1252	RotateInstrumentComponent(wksp+"det","DetectorBench",X="-1.0",Angle=str(2.0*float(incAngles[k])))
1253	Divide(LHSWorkspace=wksp+"det",RHSWorkspace=wksp+"mon",OutputWorkspace=wksp+"norm")
1254	if dlist[k] != "none":
1255	Divide(LHSWorkspace=wksp+"norm",RHSWorkspace=dlist[k],OutputWorkspace=wksp+"norm")
1256	ReplaceSpecialValues(wksp+"norm","0.0","0.0","0.0","0.0",OutputWorkspace=wksp+"norm")
1257	ConvertUnits(wksp+"norm",Target="MomentumTransfer",OutputWorkspace=wksp+"RvQ")
1258	if(diagnostics==0):
1259	DeleteWorkspace(wksp+"mon")
1260	DeleteWorkspace(wksp+"det")
1261	else:
1262	CropWorkspace(InputWorkspace=wksp,OutputWorkspace=wksp+"det",StartWorkspaceIndex=4,EndWorkspaceIndex=243)
1263	# move the first spectrum in the list onto the beam centre so that when the bench is rotated it's in the right place
1264	MoveInstrumentComponent(wksp+"det","DetectorBench",Y=str((125.0-float(minSpec))*1.2e-3))
1265	# add a bit to the angle to put the first spectrum of the group in the right place
1266	a1=2.0float(incAngles[k])+atan((float(minSpec)-float(specChan))1.2e-3/3.53)*180.0/pi
1267	#print str(2.0float(incAngles[k]))+" "+str(atan((float(minSpec)-float(specChan))1.2e-3/3.63)*180.0/pi)+" "+str(a1)
1268	RotateInstrumentComponent(wksp+"det","DetectorBench",X="-1.0",Angle=str(a1))
1269	floodnorm(wksp+"det",floodfile)
1270	if (subbgd==1):
1271	# Subract a per spectrum background
1272	Minus(wksp+"det",wbgdtemp,OutputWorkspace=wksp+"det")
1273	ResetNegatives(InputWorkspace=wksp+"det",OutputWorkspace=wksp+"det",AddMinimum='0',ResetValue="0.0")
1274	GroupDetectors(wksp+"det",OutputWorkspace=wksp+"sum",WorkspaceIndexList=range(int(minSpec)-5,int(maxSpec)-5+1),KeepUngroupedSpectra="0")
1275	else:
1276	GroupDetectors(wksp+"det",OutputWorkspace=wksp+"sum",WorkspaceIndexList=range(int(minSpec)-5,int(maxSpec)-5+1),KeepUngroupedSpectra="0")
1277	RebinToWorkspace(WorkspaceToRebin=wksp+"sum",WorkspaceToMatch=wksp+"mon",OutputWorkspace=wksp+"sum")
1278	Divide(LHSWorkspace=wksp+"sum",RHSWorkspace=wksp+"mon",OutputWorkspace=wksp+"norm")
1279	RebinToWorkspace(WorkspaceToRebin=wksp+"det",WorkspaceToMatch=wksp+"mon",OutputWorkspace=wksp+"det")
1280	Divide(LHSWorkspace=wksp+"det",RHSWorkspace=wksp+"mon",OutputWorkspace=wksp+"detnorm")
1281	if dlist[k] != "none":
1282	Divide(LHSWorkspace=wksp+"norm",RHSWorkspace=dlist[k],OutputWorkspace=wksp+"norm")
1283	ReplaceSpecialValues(wksp+"norm","0.0","0.0","0.0","0.0",OutputWorkspace=wksp+"norm")
1284	Divide(LHSWorkspace=wksp+"detnorm",RHSWorkspace=dlist[k],OutputWorkspace=wksp+"detnorm")
1285	ReplaceSpecialValues(wksp+"detnorm","0.0","0.0","0.0","0.0",OutputWorkspace=wksp+"detnorm")
1286	ConvertUnits(wksp+"norm",OutputWorkspace=wksp+"RvQ",Target="MomentumTransfer")
1287	if(diagnostics == 0):
1288	DeleteWorkspace(wksp+"sum")
1289	DeleteWorkspace(wksp+"mon")
1290	DeleteWorkspace(wksp+"det")
1291	if doCorrs != "0":
1292	if nper == 2:
1293	nrPNRCorrection(i+"norm_"+pnums[0],i+"norm_"+pnums[1],calibration=calibration)
1294	for j in range(2):
1295	RenameWorkspace(InputWorkspace=i+"norm"+"_"+pnums[j]+"corr",OutputWorkspace=i+"normcorr"+"_"+pnums[j])
1296	GroupWorkspaces(InputWorkspaces=i+"normcorr_"+pnums[0]+","+i+"normcorr_"+pnums[1],OutputWorkspace=i+"normcorr")
1297	ConvertUnits(InputWorkspace=i+"normcorr",OutputWorkspace=i+"normcorrRvQ",Target="MomentumTransfer")
1298	if (nspec > 4 and doLDCorrs != "0"):
1299	nrPNRCorrection(i+"detnorm_"+pnums[0],i+"detnorm_"+pnums[1],calibration=calibration)
1300	for j in range(4):
1301	RenameWorkspace(InputWorkspace=i+"detnorm"+"_"+pnums[j]+"corr",OutputWorkspace=i+"detnormcorr"+"_"+pnums[j])
1302	GroupWorkspaces(InputWorkspaces=i+"detnormcorr_"+pnums[0]+","+i+"detnormcorr_"+pnums[1],OutputWorkspace=i+"detnormcorr")
1303	else:
1304	nrPACorrection(i+"norm"+"_"+pnums[0],i+"norm"+"_"+pnums[1],i+"norm"+"_"+pnums[2],i+"norm"+"_"+pnums[3],calibration=calibration)
1305	for j in range(4):
1306	RenameWorkspace(InputWorkspace=i+"norm"+"_"+pnums[j]+"corr",OutputWorkspace=i+"normcorr"+"_"+pnums[j])
1307	GroupWorkspaces(InputWorkspaces=i+"normcorr_"+pnums[0]+","+i+"normcorr_"+pnums[1]+","+i+"normcorr_"+pnums[2]+","+i+"normcorr_"+pnums[3]+"",OutputWorkspace=i+"normcorr")
1308	ConvertUnits(InputWorkspace=i+"normcorr",OutputWorkspace=i+"normcorrRvQ",Target="MomentumTransfer")
1309	if (nspec > 4 and doLDCorrs != "0"):
1310	nrPACorrection(i+"detnorm_"+pnums[0],i+"detnorm_"+pnums[1],i+"detnorm_"+pnums[2],i+"detnorm_"+pnums[3],calibration=calibration)
1311	for j in range(4):
1312	RenameWorkspace(InputWorkspace=i+"detnorm"+"_"+pnums[j]+"corr",OutputWorkspace=i+"detnormcorr"+"_"+pnums[j])
1313	GroupWorkspaces(InputWorkspaces=i+"detnormcorr_"+pnums[0]+","+i+"detnormcorr_"+pnums[1]+","+i+"detnormcorr_"+pnums[2]+","+i+"detnormcorr_"+pnums[3]+"",OutputWorkspace=i+"detnormcorr")
1314	if (diagnostics == 0 and doCorrs != "0"):
1315	DeleteWorkspace(i+"norm")
1316	#DeleteWorkspace(i+"RvQ")
1317	if (diagnostics == 0 and doLDCorrs != "0"):
1318	DeleteWorkspace(i+"detnorm")
1319	k=k+1
1320	DeleteWorkspace(i)
1321	if (subbgd==1):
1322	DeleteWorkspace("wbgdtemp")
1323	#
1324	#===========================================================
1325	#
1326	def tl(wksp,th0,schan):
1327	pixel=1.2
1328	dist=3630
1329	ThetaInc=th0*pi/180.0
1330	a1=mtd[wksp]
1331	y=a1.readY(0)
1332	ntc=len(y)
1333	nspec=a1.getNumberHistograms()
1334	x1=n.zeros((nspec+1,ntc+1))
1335	theta=n.zeros((nspec+1,ntc+1))
1336	y1=n.zeros((nspec,ntc))
1337	e1=n.zeros((nspec,ntc))
1338	for i in range(0,nspec):
1339	x=a1.readX(i)
1340	y=a1.readY(i)
1341	e=a1.readE(i)
1342	x1[i,0:ntc+1]=x[0:ntc+1]
1343	theta[i,:]=atan2( (i - schan-0.5) * pixel + dist * tan(ThetaInc) , dist)*180/pi
1344	y1[i,0:ntc]=y[0:ntc]
1345	e1[i,0:ntc]=e[0:ntc]
1346	x1[nspec,:]=x1[nspec-1,:]
1347	theta[nspec,:]=atan2( (nspec - schan-0.5) * pixel + dist * tan(ThetaInc) , dist)*180/pi
1348	d1=[x1,theta,y1,e1]
1349	return d1
1350	#
1351	#===========================================================
1352	#
1353	def writemap_tab(dat,th0,spchan,fname):
1354	a1=tl(dat,th0,spchan)
1355	f=open(fname,'w')
1356	x=a1[0]
1357	y=a1[1]
1358	z=a1[2]
1359	e=a1[3]
1360	s="\t"
1361	for i in range(0,n.shape(z)[1]-1):
1362	s+="%g\t" % ((x[0][i]+x[0][i+1])/2.0)
1363	s+="%g\t" % ((x[0][i]+x[0][i+1])/2.0)
1364	s+="\n"
1365	f.write(s)
1366	for i in range(0,n.shape(y)[0]-1):
1367	s=""
1368	s+="%g\t" % ((y[i][0]+y[i+1][0])/2.0)
1369	for j in range(0,n.shape(z)[1]-1):
1370	s+="%g\t" % z[i][j]
1371	s+="%g\t" % e[i][j]
1372	s+="\n"
1373	f.write(s)
1374	f.close()
1375	#
1376	#===========================================================
1377	#
1378	def xye(wksp):
1379	a1=mtd[wksp]
1380	x1=a1.readX(0)
1381	X1=n.zeros((len(x1)-1))
1382	for i in range(0,len(x1)-1):
1383	X1[i]=(x1[i]+x1[i+1])/2.0
1384	y1=a1.readY(0)
1385	e1=a1.readE(0)
1386	d1=[X1,y1,e1]
1387	return d1
1388	#
1389	#======================================================================================
1390	#
1391	def writeXYE_tab(dat,fname):
1392	a1=xye(dat)
1393	f=open(fname,'w')
1394	x=a1[0]
1395	y=a1[1]
1396	e=a1[2]
1397	s=""
1398	s+="x\ty\te\n"
1399	f.write(s)
1400	for i in range(len(x)):
1401	s=""
1402	s+="%f\t" % x[i]
1403	s+="%f\t" % y[i]
1404	s+="%f\n" % e[i]
1405	f.write(s)
1406	f.close()
1407
1408
1409	'''
1410	def quickPlot(runlist,dataDir,lmin,reb,lmax,spmin,spmax,output,plotper,polper,zmin,zmax,zlog):
1411
1412	isisDataDir=dataDir
1413	logger.notice("setting dataDir="+dataDir+" "+isisDataDir)
1414	deleteFromRootName(output)
1415	addruns(runlist,output)
1416	nper=nperFromList(output)
1417	rebpars=str(lmin)+","+str(reb)+","+str(lmax)
1418
1419	if nper == 1:
1420	ConvertUnits(InputWorkspace=output,OutputWorkspace=output,Target="Wavelength",AlignBins="1")
1421	Rebin(InputWorkspace=output,OutputWorkspace=output,Params=rebpars)
1422	CropWorkspace(InputWorkspace=output,OutputWorkspace=output+"m",StartWorkspaceIndex="1",EndWorkspaceIndex="1")
1423	CropWorkspace(InputWorkspace=output,OutputWorkspace=output+"d",StartWorkspaceIndex=str(spmin+4),EndWorkspaceIndex=str(spmax+4))
1424	Divide(output+"d",output+"m",OutputWorkspace=output+"n")
1425	workspace_mtx=mantidplot.importMatrixWorkspace(output+"n")
1426	gr2d=workspace_mtx.plotGraph2D()
1427	l=gr2d.activeLayer()
1428	if zlog == 0:
1429	l.setAxisScale(1,zmin,zmax,0)
1430	else:
1431	l.setAxisScale(1,zmin,zmax,1)
1432	else:
1433	workspace_mtx=[]
1434	nplot=0
1435	for i in plotper:
1436	ConvertUnits(InputWorkspace=output+"_"+str(i),OutputWorkspace=output+"_"+str(i),Target="Wavelength",AlignBins="1")
1437	Rebin(InputWorkspace=output+"_"+str(i),OutputWorkspace=output+"_"+str(i),Params=rebpars)
1438	CropWorkspace(InputWorkspace=output+"_"+str(i),OutputWorkspace=output+"_"+str(i)+"m",StartWorkspaceIndex="1",EndWorkspaceIndex="1")
1439	CropWorkspace(InputWorkspace=output+"_"+str(i),OutputWorkspace=output+"_"+str(i)+"d",StartWorkspaceIndex=str(spmin+4),EndWorkspaceIndex=str(spmax+4))
1440	Divide(output+"_"+str(i)+"d",output+"_"+str(i)+"m",OutputWorkspace=output+"_"+str(i)+"n")
1441	workspace_mtx.append(mantidplot.importMatrixWorkspace(output+"_"+str(i)+"n"))
1442	gr2d=workspace_mtx[nplot].plotGraph2D()
1443	nplot=nplot+1
1444	l=gr2d.activeLayer()
1445	if zlog == 0:
1446	l.setAxisScale(1,zmin,zmax,0)
1447	else:
1448	l.setAxisScale(1,zmin,zmax,1)
1449
1450	up=mtd[output+"_2d"]
1451	down=mtd[output+"_1d"]
1452	asym=(up-down)/(down+up)
1453	RenameWorkspace(asym,OutputWorkspace=output+"_asym")
1454	ReplaceSpecialValues(output+"_asym",OutputWorkspace=output+"_asym","0.0","0.0","0.0","0.0")
1455	workspace_mtx.append(mantidplot.importMatrixWorkspace(output+"_asym"))
1456	gr2d=workspace_mtx[nplot].plotGraph2D()
1457	l=gr2d.activeLayer()
1458	l.setAxisScale(1,-1.0,1.0,0)
1459
1460
1461	logger.notice("quickPlot Finished)
1462	'''
1463

Download in other formats:

Original Format