| 1 | """A test script to compare a simulation of MERLIN convolved with a foreground model |
|---|
| 2 | """ |
|---|
| 3 | from mantid.simpleapi import * |
|---|
| 4 | from mantid.kernel import V3D |
|---|
| 5 | import os |
|---|
| 6 | |
|---|
| 7 | def create_cuboid_xml(xlength,ylength,zlength): |
|---|
| 8 | xml = """<cuboid id="sample0"> |
|---|
| 9 | <left-front-bottom-point x="%(xpt)f" y="-%(ypt)f" z="-%(zpt)f" /> |
|---|
| 10 | <left-front-top-point x="%(xpt)f" y="-%(ypt)f" z="%(zpt)f" /> |
|---|
| 11 | <left-back-bottom-point x="-%(xpt)f" y="-%(ypt)f" z="-%(zpt)f" /> |
|---|
| 12 | <right-front-bottom-point x="%(xpt)f" y="%(ypt)f" z="-%(zpt)f" /> |
|---|
| 13 | </cuboid> |
|---|
| 14 | <algebra val="sample0" /> |
|---|
| 15 | """ |
|---|
| 16 | return xml % {"xpt": xlength/2.0,"ypt":ylength/2.0,"zpt":zlength/2.0} |
|---|
| 17 | |
|---|
| 18 | |
|---|
| 19 | data_dir = '/home/dmn58364/mantidproject/test-scripts/inelastic/tobyfit/mantid' |
|---|
| 20 | |
|---|
| 21 | ei = 300. |
|---|
| 22 | bins = [-30,3,279] |
|---|
| 23 | temperature = 6. |
|---|
| 24 | chopper_speed = 600. |
|---|
| 25 | |
|---|
| 26 | # Oriented lattice & goniometer |
|---|
| 27 | alatt = 5.57 |
|---|
| 28 | blatt = 5.51 |
|---|
| 29 | clatt = 12.298 |
|---|
| 30 | uvec = [9.700000e-03,9.800000e-03,9.996000e-01] |
|---|
| 31 | vvec = [9.992000e-01,-3.460000e-02,-4.580000e-02] |
|---|
| 32 | psi = 0.0 |
|---|
| 33 | omega = 0.0 |
|---|
| 34 | alpha = 0.0 |
|---|
| 35 | beta = 0.0 |
|---|
| 36 | gamma = 0.0 |
|---|
| 37 | |
|---|
| 38 | # sample dimensions |
|---|
| 39 | sx = 0.05 |
|---|
| 40 | sy = 0.025 |
|---|
| 41 | sz = 0.04 |
|---|
| 42 | |
|---|
| 43 | # Crystal mosaic |
|---|
| 44 | eta_sig = 4.0 |
|---|
| 45 | |
|---|
| 46 | ## |
|---|
| 47 | ## A fake data set for simulation purposes. |
|---|
| 48 | ## |
|---|
| 49 | fake_data = CreateSimulationWorkspace(Instrument='MERLIN',BinParams=bins,UnitX='DeltaE', |
|---|
| 50 | DetectorTableFilename=os.path.join(data_dir,'MER06398.raw')) |
|---|
| 51 | |
|---|
| 52 | nhist = fake_data.getNumberHistograms() |
|---|
| 53 | nbins = fake_data.blocksize() |
|---|
| 54 | print 'nbins*nhist',nhist*nbins |
|---|
| 55 | |
|---|
| 56 | ## |
|---|
| 57 | ## Required log entries, can be taken from real ones by placing an instrument parameter of the same |
|---|
| 58 | ## name pointing to the log name |
|---|
| 59 | ## |
|---|
| 60 | AddSampleLog(Workspace=fake_data, LogName='Ei',LogText=str(ei), LogType="Number") |
|---|
| 61 | AddSampleLog(Workspace=fake_data, LogName='temperature_log',LogText=str(temperature), LogType="Number") |
|---|
| 62 | AddSampleLog(Workspace=fake_data, LogName='chopper_speed_log',LogText=str(chopper_speed), LogType="Number") |
|---|
| 63 | AddSampleLog(Workspace=fake_data, LogName='eta_sigma',LogText=str(eta_sig), LogType="Number") |
|---|
| 64 | |
|---|
| 65 | |
|---|
| 66 | ## |
|---|
| 67 | ## OrientedLattice & Goniometer to define rotations |
|---|
| 68 | ## |
|---|
| 69 | SetUB(Workspace=fake_data,a=alatt,b=blatt,c=clatt,u=uvec,v=vvec) |
|---|
| 70 | # AddSampleLog(Workspace=fake_data,LogName='psi',LogText=str(psi),LogType='Number') |
|---|
| 71 | # AddSampleLog(Workspace=fake_data,LogName='omega',LogText=str(omega),LogType='Number') |
|---|
| 72 | # AddSampleLog(Workspace=fake_data,LogName='alpha',LogText=str(alpha),LogType='Number') |
|---|
| 73 | # AddSampleLog(Workspace=fake_data,LogName='beta',LogText=str(beta),LogType='Number') |
|---|
| 74 | # AddSampleLog(Workspace=fake_data,LogName='gamma',LogText=str(gamma),LogType='Number') |
|---|
| 75 | # SetGoniometer(Workspace=fake_data,Axis0="psi,0,1,0,1",Axis1="omega,0,1,0,1", |
|---|
| 76 | # Axis2="gamma,0,1,0,1",Axis3="beta,1,0,0,1",Axis4="alpha,0,0,1,1") |
|---|
| 77 | |
|---|
| 78 | ## |
|---|
| 79 | ## Sample shape |
|---|
| 80 | ## |
|---|
| 81 | CreateSampleShape(InputWorkspace=fake_data, ShapeXML=create_cuboid_xml(sx,sy,sz)) |
|---|
| 82 | |
|---|
| 83 | ## |
|---|
| 84 | ## Chopper & Moderator models. |
|---|
| 85 | ## |
|---|
| 86 | CreateModeratorModel(Workspace=fake_data,ModelType='IkedaCarpenterModerator', |
|---|
| 87 | Parameters="TiltAngle=32,TauF=2.7,TauS=0,R=0") |
|---|
| 88 | CreateChopperModel(Workspace=fake_data,ModelType='FermiChopperModel', |
|---|
| 89 | Parameters="AngularVelocity=chopper_speed_log,ChopperRadius=0.049,SlitThickness=0.0023,SlitRadius=1.3,Ei=Ei") |
|---|
| 90 | |
|---|
| 91 | ## |
|---|
| 92 | ## Create the MD workspace |
|---|
| 93 | ## |
|---|
| 94 | simulws = ConvertToMD(InputWorkspace=fake_data, QDimensions="Q3D", SplitInto=[3], SplitThreshold=100, |
|---|
| 95 | MinValues="-15,-15,-15,-30", MaxValues="25,25,25,279") |
|---|
| 96 | |
|---|
| 97 | ## |
|---|
| 98 | ## Run the simulation |
|---|
| 99 | ## |
|---|
| 100 | resol_model = "TobyFitResolutionModel" |
|---|
| 101 | xsec_model = "Strontium122" |
|---|
| 102 | parameters = "Seff=0.7,J1a=38.7,J1b=-5.0,J2=27.3,SJc=10.0,GammaSlope=0.08,MultEps=0,TwinType=0,MCLoopMax=1000" |
|---|
| 103 | simulated_ws = SimulateResolutionConvolvedModel(InputWorkspace=simulws, |
|---|
| 104 | ResolutionFunction=resol_model, |
|---|
| 105 | ForegroundModel=xsec_model, |
|---|
| 106 | Parameters=parameters) |
|---|
| 107 | |
|---|
| 108 | #simulated_file = 'mer300meV_psi0_sim_mantid_1.nxs' |
|---|
| 109 | #SaveMD(InputWorkspace=simulated_ws,Filename=os.path.join(data_dir,simulated_file)) |
|---|
