C-- SANE_HMS_physics from HMS data C-- ebeam = (GeV) Energy of incident beam C-- HMS_mom = hsdelta(GeV) Central momentum of HMS C-- HMS_angle = hstheta(radian) Angle between HMS and outgoing beam c Written to work with hms[run].rzdat files real function SANE_HMS_pf() implicit none c include 'sane.cmn' include ? logical first data first /.true./ save first real ebeam,bigcmom,bigcang c real ebigc,thbigc c real xbSANE,Q2SANE,W2SANE,WSANE c real nuSANE,ESANE,thetaSANE,ySANE real xbHMS,Q2HMS,W2HMS,WHMS real nuHMS,EHMS,dEHMS,thetaHMS,phiHMS,yHMS real targxp,targyp,targx,targy,beamz,SRy,SRx real targxpfp,targypfp,targxfp,targyfp real mpr,d2r,r2d,hel real wt,sc_fp cc real trkeff(16) cc data trkeff/1.000,0.994,0.998,0.998,0.997,0.998,0.999,0.999, cc & 0.996,0.998,0.999,0.999,0.997,0.999,1.000,1.000/ integer i,si integer n_W real Wmin,Wmax,Emin,Emax vector trkeff vector trkefffit vector runinfo vector Winfo mpr = 0.938272 d2r = 3.14159/180. r2d = 1 / d2r ccc ebeam = 4.734 cc ebeam = 4.731 c ebeam = 5.892 cc ebeam = 5.893 cc ebeam = 5.896 ebeam = runinfo(1)/1000. n_W=int(Winfo(1)) Wmin=Winfo(2) Wmax=Winfo(3) Emin=runinfo(2)-0.4 Emax=runinfo(2)+0.4 if(first) then first=.false. call hbook2(100,'hsxfp vs. W',50,0.,4.0,50,-50.,50.,0) call hbook1(10000,'Q^2! ' ,50,0.5,4.,0) call hbook1(10001,'Q^2! ' ,50,0.,7.,0) call hbook1(20000,'x?bj! ' ,50,0.,1.,0) call hbook1(20001,'x?bj! ' ,50,0.,2.,0) c call hbook1(300000,'W ' ,100,0.7,3.0,0) cc call hbook1(30001,'W ' ,50,0.7,1.8,0) call hbook1(30000,'W ' ,21,2.1475,2.2525,0) call hbook1(30000,'W ' ,27,1.8475,1.9825,0) call hbook1(30000,'W ' ,81,1.1975,1.6025,0) call hbook1(30000,'W ' ,100,0.9,2.9,0) cc call hbook1(30000,'W ' ,50,0.9,1.8,0)!CHe-C cc call hbook1(30000,'W ' ,35,1.090,1.615,0)!4.7_p3200th20 call hbook1(30000,'W ',n_W,Wmin,Wmax,0) call hbook1(40000,'E?HMS! ' ,50,Emin,Emax,0) cc call hbook1(30000,'W ' ,13,2.1,2.25,0)!CHe-C cc call hbook1(30000,'W ' ,50,1.9,2.5,0)!CHe-C cc call hbook1(40000,'E?HMS! ' ,50,1.7,2.5,0) cc call hbook1(40000,'E?HMS! ' ,50,2.2,3.,0) cc call hbook1(40000,'E?HMS! ' ,50,2.7,3.5,0) cc call hbook1(40000,'E?HMS! ' ,50,2.5,4.,0) c call hbook1(40000,'E?HMS! ' ,50,3.8,5.,0) call hbook1(41000,'[D]E?HMS! ' ,50,-18.,18.,0) c call hbook1(50000,'theta?HMS! ' ,50,0.,40.,0) call hbook1(50000,'[Q]?HMS! ' ,50,12.,25.,0) c call hbook1(50000,'[Q]?HMS! ' ,50,15.,25.,0) call hbook1(60000,'nu?HMS! ' ,50,0.,4.,0) call hbook1(70000,'y ' ,50,0.,1.,0) call hbook2(120000,'Q^2! vs. x?bj! HMS' ,50,0.,2., > 50,0.,7.0,0) c call hbook2(130000,'Q^2! vs. W HMS' ,50,0.1,4., c > 50,0.,7.0,0) call hbook2(130000,'Q^2! vs. W HMS' ,50,0.5,2., > 50,0.,4.0,0) call hbook2(480000,'E vs. [Q] HMS' ,50,15.,25., > 50,2.,4.0,0) call hbook1(4100,'SR?x!' ,50,-1.5,1.5,0) call hbook1(4200,'SR?y!' ,50,-1.5,1.5,0) call hbook1(1000,'xptarg' ,50,-0.23,0.13,0) call hbook1(1001,'xptarg no [D]P cut' ,50,-0.13,0.13,0) call hbook1(1100,'xptargfp' ,50,-0.08,0.08,0) call hbook1(1200,'xtargfp' ,50,-50.,50.,0) call hbook1(2000,'yptarg' ,50,-0.06,0.06,0) call hbook1(2001,'yptarg no [D]P cut' ,50,-0.06,0.06,0) call hbook1(2100,'yptargfp' ,50,-0.04,0.03,0) call hbook1(2200,'ytargfp' ,50,-30.,25.,0) call hbook1(3000,'beam z' ,50,-10.,10.,0) call hbook1(3001,'beam z no [D]P cut' ,50,-10.,10.,0) call hbook2(12000,'xptar vs. yptar' ,50,-0.06,0.06, > 50,-0.23,0.13,0) call hbook2(141000,'xptar vs. [D]E?HMS!' ,50,-9.,9., > 50,-0.23,0.13,0) call hbook2(111100,'xptar vs. xtar' ,50,-1.5,1.5, > 50,-0.13,0.13,0) call hbook2(211100,'ytar vs. xtar' ,50,-1.5,1.5, > 50,-4.,4.,0) call hbook2(311000,'beam z vs. xtar' ,50,-1.5,1.5, > 50,-10.,10.,0) call hbook2(321000,'beam z vs. ytar' ,50,-4.,4., > 50,-10.,10.,0) call hbook2(412000,'SR?x! vs. ytar' ,50,-2.5,2.5, > 50,-1.5,1.5,0) call hbook2(412100,'SR?x! vs. yptar' ,50,-0.06,0.06, > 50,-1.5,1.5,0) call hbook2(421100,'SR?y! vs. xptar' ,50,-0.23,0.13, > 50,-1.5,1.5,0) call hbook1(1110,'xtar',50,-4.,4.,0) call hbook1(1111,'ytar',50,-4.,4.,0) call hbook1(1112,'[F]?HMS! ' ,100,120,200.,0) call hbook2(1113,'[D]E?HMS vs W',n_W,Wmin,Wmax,50,-18.,18.,0) endif C do i=1,nclust C C hel = i_helicity C ebigc = eclust(i)*1000 C nuSANE = ebeam - ebigc C thetaSANE = r2d*atan(xclust(i)/345.) + 40. C Q2SANE = 4.* ebeam*ebigc * sin(d2r*thetaSANE/2)**2 C xbSANE = Q2SANE / 2. / mpr / nuSANE C W2SANE = mpr**2 + 2.*mpr*nuSANE - Q2SANE C WSANE = sqrt(max(0.,W2SANE)) C ySANE = 1 - ebigc/ebeam C C enddo dEHMS = hsdelta EHMS = hse c EHMS = hse-0.0016 c EHMS = hse-0.0037 This seems to flatten out W dat/MC ratio for 73027, w/o shifting theta c EHMS = hse-0.0016 For fitting W from 1.2-1.6GeV (Delta Res.) nuHMS = ebeam - EHMS thetaHMS = hstheta*r2d thetaHMS = thetaHMS!-0.018 phiHMS = hsphi*r2d c thetaHMS = thetaHMS-0.055 This seems to flatten out W dat/MC ratio for 73027, w/o shifting E' c thetaHMS = thetaHMS-0.018 For fitting W from 1.2-1.6GeV (Delta Res.) Para c thetaHMS = thetaHMS-0.075 For fitting W from 1.2-1.6GeV (Delta Res.) Perp c thetaHMS = thetaHMS-0.60 W from 1.85-1.98GeV Perp c thetaHMS = thetaHMS-0.095 5.9 Perp 4.4GeV/c 15.4degrees c Q2HMS = 4.* ebeam*EHMS * sin(d2r*thetaHMS/2)**2 c Q2HMS = 4.* ebeam*EHMS * sin(hstheta/2)**2 Q2HMS = 4.* ebeam*EHMS * sin((thetaHMS*d2r)/2)**2 xbHMS = Q2HMS / 2. / mpr / nuHMS W2HMS = mpr**2 + 2.*mpr*nuHMS - Q2HMS WHMS = sqrt(max(0.,W2HMS)) c WHMS = WHMS+0.00765 !!! w/o beam offset w/ overall z shift of 1.5mm + 2.5mm for C c WHMS = WHMS+0.0051 !!! W shift w/ beampos/SR offsets no dE on MC 06/11/09 c WHMS = WHMS+0.012 !!! W shift w/o beampos/SR offsets 06/11/09 c WHMS = WHMS+0.0113 !!! W shift with beampos/SR offsets 06/10/09 c WHMS = WHMS+0.0021 !!! W shift with SR correction to hszbeam 05/28/09 c WHMS = WHMS+0.0035 !!! W shift without SR correction to hszbeam 05/27/09 yHMS = 1 - EHMS/ebeam targxp = hsxptar targyp = hsyptar targx = hsxtar targy = hsytar beamz = hszbeam SRy = srast_y SRx = srast_x targxpfp = hsxpfp targypfp = hsypfp targxfp = hsxfp targyfp = hsyfp c print*,Q2 wt=1. cc sc_fp=33.750 cc do si=1,10 cc if(abs(hsxfp-sc_fp).le.3.75) then cc wt=1./trkeff(si) cc write(*,*) hsxfp, sc_fp, trkeff(si) cc endif cc sc_fp=sc_fp-7.5 cc enddo ccc if(abs(hsxfp).lt.100) then ccc wt=1./(trkefffit(3)*hsxfp**2+trkefffit(2)*hsxfp ccc & +trkefffit(1)) ccc endif CCC if(abs(hsxfp).lt.100) then CCC wt=1./(trkefffit(4)*hsxfp**3+!trkefffit(5)*hsxfp**4+trkefffit(4)*hsxfp**3+ CCC & trkefffit(3)*hsxfp**2+trkefffit(2)*hsxfp+trkefffit(1)) CCC endif ccc if(wt.gt.1./0.5) then ccc wt=0. ccc endif CCC You need to turn off this part if you want to use fitted one ccc if(num_1xsc.ge.4.and.num_1xsc.le.13) then ccc wt=1./trkeff(num_1xsc-3) ! Inintial scin # = 4 ccc else ccc wt=0. ccc endif CCC if(num_1xsc.eq.9)then CCC wt=1.2 CCC endif ccc wt=1. if(WHMS.ge.Wmin.and.WHMS.le.Wmax)then if(hse.gt.0.and.(hsshtrk/hse).gt.0.7)then c if(abs(hsxfp).le.60)then c if(thetaHMS.gt.16.and.thetaHMS.lt.20.and.abs(dEHMS).lt.10) then ccc if(dEHMS.lt.-1.and.dEHMS.gt.-30)then c if(abs(dEHMS).lt.30.and.abs(hsxfp).lt.37.5)then !.and.abs(dEHMS).gt.1.!.and.abs(hsxfp).lt.30 CCC if(num_1xsc.ge.4.and.num_1xsc.le.13)then ccc if(dEHMS.gt.-8.and.dEHMS.lt.12)then if(dEHMS.gt.-8.and.dEHMS.lt.8)then ccc if(EHMS.lt.runinfo(5)*1.08.and.EHMS.gt.runinfo(5)*0.92)then ccc & .and.SRy.gt.0.)then call hf1(3001,beamz,wt) call hf1(41000,dEHMS,wt) call hf2(141000,dEHMS,targxp,wt) c if(targxp.gt.-0.2.and.targxp.lt.-0.1)then cc if(cnt_hk.ge.1.and.cnt_hk.le.16) then cc wt=1./trkeff(cnt_hk) cc else cc wt=1. cc endif call hf1(1000,targxp,wt) call hf1(2000,targyp,wt) call hf1(1100,targxpfp,wt) call hf1(2100,targypfp,wt) call hf1(1200,targxfp,wt) call hf1(2200,targyfp,wt) call hf1(3000,beamz,wt) call hf1(4100,SRx,wt) call hf1(4200,SRy,wt) call hf2(412000,targy,SRx,wt) call hf2(412100,targyp,SRx,wt) call hf2(421100,targxp,SRy,wt) call hf1(10000,Q2HMS,wt) call hf1(20000,xbHMS,wt) call hf1(300000,WHMS,wt) call hf1(30000,WHMS,wt) call hf1(40000,EHMS,wt) call hf1(50000,thetaHMS,wt) call hf1(60000,nuHMS,wt) call hf1(70000,yHMS,wt) call hf2(12000,targyp,targxp,wt) call hf2(111100,targx,targxp,wt) call hf2(211100,targx,targy,wt) call hf2(311000,targx,beamz,wt) call hf2(321000,targy,beamz,wt) call hf2(120000,xbHMS,Q2HMS,wt) call hf2(130000,WHMS,Q2HMS,wt) call hf2(480000,thetaHMS,EHMS,wt) call hf2(100,WHMS,hsxfp,1.) call hf1(1110,targx,wt) call hf1(1111,targy,wt) call hf1(1112,phiHMS,wt) call hf2(1113,WHMS,dEHMS,wt) endif endif endif c endif c endif end