macro norm_correct sys=8 part=pim rapbin=0

vec/del *
hi/del *

if ([rapbin].lt.0) then
  ybin=b$eval(-1*[rapbin])
elseif ([rapbin].eq.0) then
  ybin=y[rapbin]
elseif ([rapbin].gt.0) then
  ybin=f[rapbin]
endif

vec/cre [part](2,40) r
vec/read [part] [sys]gev_[part]_[ybin]_raw.dat
vec/cre amp(40) r
vec/cre ampERR(40) r
vec/cre mt-m0(40) r

do i=1,40
vec/copy [part](1,[i]) amp([i])
vec/copy [part](2,[i]) ampERR([i])
vec/in mt-m0([i]) $eval(([i]*0.025)-0.0125)
enddo

1d 2000 ''//[sys]//'GeV '//[part]//' with new kaons' 40 0.0 1.0

*************************************************
* Normalize the results by:                     *
*                                               *
* 40 bins/unit of mt-m0                         *
* 10 bins/unit of y                             *
* 1/2pi                                         *
* 1/mt of the bin                               *
* per Event 2 GeV: 19054                        *
*           4 GeV: 21502                        *
*           6 GeV: 28385                        *
*           8 GeV: 17918                        *
*                                               *
*************************************************

* Compute mt for each mt-m0 bin:
vec/cre mt(40) r
vec/op/vbias mt-m0 0.13956995 mt

if ([sys].eq.2) then
  evtscale=0.0033411345
elseif ([sys].eq.4) then
  evtscale=0.0029607468
elseif ([sys].eq.6) then
  evtscale=0.0022428035
elseif ([sys].eq.8) then
  evtscale=0.0035529622
endif
* Normalize bins, 1/2pi, per Event
do i=1,40
 vec/op/vscale amp([i]) [evtscale] amp([i])
 vec/op/vscale ampERR([i]) [evtscale] ampERR([i])
enddo
* Normalize 1/mt
do i=1,40
   vec/op/vdivide amp([i]) mt([i]) amp([i])
   vec/op/vdivide ampERR([i]) mt([i]) ampERR([i])
enddo
*Add loss correction:
vec/cre loss(31,40) r
vec/read loss ../accept/[sys]gev_pim_full_loss.dat
vec/cre spec_loss(40) r
do i=1,40
 vec/copy loss($eval(16+[rapbin]),[i]) spec_loss([i])
enddo
vec/del loss
do i=1,40
 loss=spec_loss([i])
 if ([loss].gt.0) then
   vec/op/vdivide amp([i]) spec_loss([i]) amp([i])
   vec/op/vdivide ampERR([i]) spec_loss([i]) ampERR([i])
 else
   vec/op/vmult amp([i]) spec_loss([i]) amp([i])
   vec/op/vmult ampERR([i]) spec_loss([i]) ampERR([i])
 endif
enddo

vec/write mt-m0,amp,ampERR [sys]gev_[part]_[ybin]_final.dat '(40(F10.4,2F10.2))' 'oc'

return