program proxi ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc c programma GEANT3 per simulazione RICH proximity focusing c l'apparato e' composto da: c - radiatore in freon c - finestra in quarzo c - proximity gap (metano) c - pad detector C C Author: R. IOMMI C Modified by: E. Cisbani C Modified by: L. Pappalardo ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc #undef CERNLIB_GEANT321_GCUNIT_INC C #include "geant321/gcunit.inc" #undef CERNLIB_GEANT321_GCONST_INC C #include "geant321/gconst.inc" external fric c parameter (ng=69024,nh=1000000) parameter (ng=1000000,nh=10000000) real XPROXI(3), XFREON(3), XQUARZO(3), XGAP(3), XDET(3), XDRIFT(3),XCERE(3) real R1FREO,R2FREO,XOFREO,ANFREO,R1DETE,R2DETE,XODETE,A1DETE,A2DETE real ROTRIC(6) real pad_dx, pad_dy,pad_gx,pad_gy,pad_sx,pad_sy real etimli integer typpa integer ifine integer rad_ty real pad_nx, pad_ny c real thetamud, phimud, pmud common/radtyp/rad_ty common/padgeo/pad_dx,pad_dy,pad_gx,pad_gy,pad_sx,pad_sy,pad_nx,pad_ny common/ancora/ifine common/pawc/h(nh) common/gcbank/q(ng) common/inpdat/plabmi,dz,sigfax,sigfay,sigfte,sigfph,plabma,typpa,plabmo common/ener/mmmr common/vin1/cf1,cf2,cfz,imud common/matsiz/XPROXI,XFREON,XQUARZO,XGAP,XDET,XDRIFT,ROTRIC,XCERE common/clas12/R1FREO,R2FREO,XOFREO,ANFREO,R1DETE,R2DETE,XODETE,A1DETE,A2DETE common/grafic/inter,kwtip common/view/teta,fi,psi,u0,v0,su,sv common/color/ICOLG common/spec/ispec common/qetab/eqe(14),qe(14) integer iranse common/riccar/nrpho,ntpho,nnttph,w,iverbo parameter (maxgamma=2000) c*** Viene inclusa anche la ricostruzione dell angolo Cherenkov common/trkpar/ivent,ipar,xp(3),x0(2),ngtot,nglost,ngamma,xg(maxgamma),yg(maxgamma), + thetag(maxgamma),egamma(maxgamma),z0g(maxgamma), + theta0g(maxgamma),nfreon,xpar,ypar,thetapar,phipar,xc,yc, + thetac(maxgamma), + phic(maxgamma),x(maxgamma),y(maxgamma),b(maxgamma), + etha(maxgamma) parameter (maxpad=1000) common/padpar/npad,ipad(maxpad),jpad(maxpad),xgp(maxpad),ygp(maxpad),ngpad(maxpad),rthpad(maxpad) typpa=8 ifine=0 c*** lettura dimensioni e parametri write (*,*) 'Read configuration file: proxi.dat' open (unit=10,file='proxi.dat',status='old') read (10,2323) c cf1 = polar angle (Deg) of the normal to the rich entrance window (clas12) c cf2 = azimuthal angle (Deg) .. (clas12) c cfz = distance origin - rich entrance window in clas12 read (10,*) cf1,cf2,cfz,imud read (10,2323) read (10,*) sigfax,sigfay read (10,2323) read (10,*) sigfte,sigfph read (10,2323) read (10,*) plabmi,plabma read (10,2323) read (10,*) XPROXI(1),XPROXI(2),XPROXI(3) read (10,2323) read (10,*) XFREON(1),XFREON(2),XFREON(3),rad_ty read (10,2323) read (10,*) R1FREO,R2FREO,XOFREO,ANFREO read (10,2323) read (10,*) XQUARZO(1),XQUARZO(2),XQUARZO(3) read (10,2323) read (10,*) XGAP(1),XGAP(2),XGAP(3) read (10,2323) read (10,*) XDET(1),XDET(2),XDET(3),XDEFF read (10,2323) read (10,*) R1DETE,R2DETE,XODETE,A1DETE,A2DETE read (10,2323) read (10,*) XDRIFT(1),XDRIFT(2),XDRIFT(3) read (10,2323) read (10,*) XCERE(1),XCERE(2),XCERE(3) read (10,2323) read (10,*) ROTRIC(1),ROTRIC(2),ROTRIC(3),ROTRIC(4),ROTRIC(5),ROTRIC(6) read (10,2323) read (10,*) pad_dx,pad_dy,pad_gx,pad_gy,pad_sx,pad_sy read (10,2323) read (10,*) pad_nx, pad_ny read (10,2323) read (10,*) inter,kwtip read (10,2323) read (10,*) teta,fi,psi,u0,v0,su,sv read (10,2323) read (10,*) ispec read (10,2323) read (10,*) ICOLG read (10,2323) read (10,*) dz,mmmr read (10,2323) read (10,*) iverbo read (10,2323) read (10,*) iranse close (unit=10) c*** clas12 init radiator position parameters (see procla.fpp) call inicod(cfz, cf1/180.*3.1415, cf2/180.*3.1415) c*** momento centrale plabmo=(plabmi+plabma)/2. c*** lunghezza apparato w=2.*(XFREON(3)+XQUARZO(3)+XGAP(3)+XDET(3)) c*** lettura tavole quantum efficiency c file originale di Evaristo: qe_original_evaristo.dat c file attualmente in uso: qe.dat ottenuto da una quantum c efficiency gaussiana nell intervallo 200-600 nm, simile c ad andamenti riportati in tesi evaristo (fig. 4.15) write (*,*) 'Read Quantum Efficiency data: qe.dat' open (unit=10,file='input/qe.dat',status='old') read(10,*) (eqe(i),qe(i), i=14,1,-1) close (unit=10) do i=1,14 qe(i)=xdeff*qe(i) enddo c*** apertura file spazio delle fasi c open (unit=150,name='s2.mud',type='old') c--------- lettura nuovo file da J. Lerose ------------------- write (*,*) 'Read Trajectiory Phase Space by Lerose: sept_1.mud (unit=15)' c open (unit=15,file='input/rich.mud',status='old') open (unit=15,file='script/rich.mud',status='old') c open (unit=16,file='script/output.mud',status='unknown') write (*,*) 'Fatto ' c------------------------------------------------------------- c*** set high time limits call timest(999999.9) call timel(etimli) write(*,*) 'Initial Time Limit (s): ',etimli c*** allocazione di memoria per ZEBRA e HBOOK write (*,*) 'Zebra and HBook Initialization' call gzebra(ng) call hlimit(-nh) * if (inter .eq. 1) then call hplint(kwtip) * end if c*** Inizializzazione PAW write (*,*) 'PAW initialization: create proxie.rz' call hropen(44,'rich','proxie.rz','N',1024,ISTAT) call hbnt(10,'proxi','D') call hbname(10,'phot',ivent, + 'ivent,ipar,px,py,pz,x0,y0,ngtot,nglost,ngamma[0,2000],xg(ngamma),yg(ngamma),'// + 'thetag(ngamma),egamma(ngamma),'// + 'z0g(ngamma),'// + 'theta0g(ngamma),nfreon,xpar,ypar,thetapar,phipar,xc,yc,thetac(ngamma),'// + 'phic(ngamma),x(ngamma),y(ngamma),b(ngamma),etha(ngamma)') ngamma=0 ivent=0 call hbname(10,'digi',npad,'npad[0,1000],ipad(npad),jpad(npad),'// + 'xgp(npad),ygp(npad),ngpad(npad),'// + 'rthpad(npad)') c call hbname(10,'mud',thetamud,phimud,pmud) c*** Inizializzazione random seed (ranlux) c call rmarut(iiseed,iin1,iin2) c iiseed=iranse c call rmarin(iiseed,iin1,iin2) c call rmarut(iiseed,iin1,iin2) write (*,*) 'Random Seed inizialization' call RLUXGO(3,iranse,0,0) c*** inizializzazione write (*,*) 'Initialize GEANT' call uginit close(90) c*** procedura sugli eventi write (*,*) 'Start the Event Loop' call grun c*** terminazione write (*,*) 'Call end of GEANT' call uglast * if (inter .eq. 1) then call igend * end if stop C 100 format (q,30a1) 2323 format (1X) end c**************************************************************************** subroutine uginit c**************************************************************************** #undef CERNLIB_GEANT321_GCTMED_INC #undef CERNLIB_GEANT321_GTTMED_INC #include "geant321/gctmed.inc" #undef CERNLIB_GEANT321_GCPHYS_INC #undef CERNLIB_GEANT321_GTPHYS_INC #include "geant321/gcphys.inc" common/spec/ispec common/view/teta,fi,psi,u0,v0,su,sv integer typpa common/grafic/inter,kwtip common/inpdat/plabmi,dz,sigfax,sigfay,sigfte,sigfph,plabma,typpa,plabmo c write(6,*) 'SUBROUTINE UGINIT' write(6,*) '*' write(6,*) '********* GEANT INITIALIZATION *********' write(6,*) '*' c*** inizializza le variabili di GEANT call ginit c*** open the data card file and tell the unit number to the ffread routine c the data card contains the two random seeds used by GEANT c open (unit=90,file='input/data-card.ff',status='old') call ffset('LINP',90) c*** attivazione fotoni Cherenkov ITCKOV = 1 call gffgo c*** inizializza DATA STRUCTURES call gzinit c*** inizializza DRAWING PACKAGE call gdinit c*** definizione e stampa materiali e particelle standard call gmate c*** descrizione apparato call proxifoc call gpart c*** sezioni d'urto e perdita di energia call gphysi end c************************************************************************** subroutine uglast c************************************************************************** #undef CERNLIB_GEANT321_GCUNIT_INC C #include "geant321/gcunit.inc" #undef CERNLIB_GEANT321_GCONST_INC C #include "geant321/gconst.inc" #undef CERNLIB_GEANT321_GCKINE_INC C #include "geant321/gckine.inc" #undef CERNLIB_GEANT321_GCFLAG_INC C #include "geant321/gcflag.inc" integer typpa real XPROXI(3),XFREON(3),XQUARZO(3),XGAP(3),XDET(3),XDRIFT(3),XCERE(3) real ROTRIC(6) common/inpdat/plabmi,dz,sigfax,sigfay,sigfte,sigfph,plabma,typpa,plabmo common/ener/mmmr common/vin1/cf1,cf2,cfz,imud common/ame/b(20) common/matsiz/XPROXI,XFREON,XQUARZO,XGAP,XDET,XDRIFT,ROTRIC,XCERE common/mong/d,IYESM,dt,ds byte nomepr(20) dimension ubf(1) close(15) close(16) c write(6,*) 'SUBROUTINE UGLAST' c*** fetch particella per output call gfpart (typpa,nomepr,itrtyp,amass,charge,tlife,ubf,nubf) c*** chiamata della routine terminale standard call glast c*** chiude rz file per paw call HLDIR('.','t') call HCDIR('//rich',' ') call hrout(0,icycle,'rich') call hrend('rich') end c*************************************************************************** subroutine gukine c*************************************************************************** #undef CERNLIB_GEANT321_GCONST_INC C #include "geant321/gconst.inc" #undef CERNLIB_GEANT321_GCKINE_INC C #include "geant321/gckine.inc" #undef CERNLIB_GEANT321_GCFLAG_INC C #include "geant321/gcflag.inc" integer typpa real XPROXI(3),XFREON(3),XQUARZO(3),XGAP(3),XDET(3),XDRIFT(3),XCERE(3) real R1FREO,R2FREO,XOFREO,ANFREO,R1DETE,R2DETE,XODETE,A1DETE,A2DETE real ROTRIC(6) real xranlu(5) real thetamud, phimud, pmud real zdummy real x_sec, y_sec, a_sec, r_sec real xf_min real ther, phir integer n_out integer ifine common/ancora/ifine common/matsiz/XPROXI,XFREON,XQUARZO,XGAP,XDET,XDRIFT,ROTRIC,XCERE common/clas12/R1FREO,R2FREO,XOFREO,ANFREO,R1DETE,R2DETE,XODETE,A1DETE,A2DETE common/inpdat/plabmi,dz,sigfax,sigfay,sigfte,sigfph,plabma,typpa,plabmo common/vin1/cf1,cf2,cfz,imud dimension vertice(3),plab(3),ub(1) real evathe, evaphi, dp SAVE plab, ub, vertice, evathe, evaphi, dp common/riccar/nrpho,ntpho,nnttph,w,iverbo parameter (maxgamma=2000) common/trkpar/ivent,ipar,xp(3),x0(2),ngtot,nglost,ngamma,xg(maxgamma),yg(maxgamma), + thetag(maxgamma),egamma(maxgamma),z0g(maxgamma), + theta0g(maxgamma),nfreon,xpar,ypar,thetapar,phipar,xc,yc, + thetac(maxgamma), + phic(maxgamma),x(maxgamma),y(maxgamma),b(maxgamma), + etha(maxgamma) REAL xin, yin, zin, xinp, yinp REAL pigre INTEGER tfla c write(6,*) 'SUBROUTINE GUKINE' sigfax = R2FREO xf_min = XOFREO sigfay = R2FREO * SIN(ANFREO/180.*3.1415) if (iverbo.gt.3) then write(6,*)' IVERBO-3 : gukine started' endif nrpho=0 ntpho=0 if(typpa.gt.14) then typpa=8 endif n_out = 0 tfla = 0 3321 continue if(typpa.eq.8) then if ((sigfte.eq.0.).and.(sigfph.eq.0.)) then c--------------- lettura nuovo file da J. Lerose --------------------- if (imud .eq. 0) then read(15,'(10e14.6)',END=1136) (xdum, i=1,4),vertice(1),evathe,vertice(2),evaphi,xdum,dp endif if (imud .eq. 1) then read(15,'(5e14.6)',END=1136) vertice(1),evathe,vertice(2),evaphi,dp endif vertice(1)=vertice(1)*100.+cf1+cfz*tan(evathe) vertice(2)=vertice(2)*100.+cf2+cfz*tan(evaphi) dp=dp+1. c clas12 if (imud .eq. 2) then read(15,'(6F10.4)',END=1136) vertice(1), vertice(2), evathe, evaphi, plabmo, zdummy c write(16,111) vertice(1),vertice(2), evathe, evaphi, plabmo, zdummy thetamud=evathe phimud=evaphi pmud=plabmo c call hbname(10,'mud',thetamud,phimud,pmud) pigre=3.14159265 dp=1. C rotate the tracks to center in sector at 0 , 60, 120, ... C mod is used to bring all track in a single sector ... increase statistics C the detector is symmetric evaphi = mod(mod(evaphi+.83,2.*pigre),pigre/3.)-pigre/6.+.2 C write(*,*) 'CLAS ref:',evathe, evaphi, xin, yin call impaxy(vertice(1), vertice(2), 0., evathe, evaphi, xin, yin, zin, ther) call xyz2fr(xin, yin, zin, xinp, yinp, phir) vertice(1) = xinp vertice(2) = yinp evathe = ther evaphi = phir tfla=1 C write(*,*) 'RICH ref:',evathe, evaphi, xinp, yinp endif c evathe=atan(evathe) c evaphi=atan(evaphi) else c--------------- random numbers call ranlux(xranlu,5) vertice(1)=sigfax*xranlu(1)+xf_min vertice(2)=2.*sigfay*(xranlu(2)-0.5) evathe=sigfte*xranlu(3) evaphi=2.*sigfph*(xranlu(4)-0.5) dp=((plabma-plabmi)*xranlu(5)+plabmi)/plabmo endif if ((iverbo .gt. 0) .and. (MOD(ivent,100) .eq. 0) .and. (ivent .GT. 0)) then write(6,*)' Track Number ',ivent/3, ' , attempts ', n_out endif c---------------------------------------------------------------------- if (tfla .EQ. 1) then plab(1)=sin(evathe)*cos(evaphi)*plabmo*dp plab(2)=sin(evathe)*sin(evaphi)*plabmo*dp plab(3)=cos(evathe)*plabmo*dp else plab(3)=plabmo*dp/sqrt(tan(evathe)*tan(evathe)+tan(evaphi)*tan(evaphi)+1.) plab(1)=tan(evathe)*plab(3) plab(2)=tan(evaphi)*plab(3) endif C neoceram ??? vertice(3)=-5. c vertice(1)=vertice(1)-dz*plab(1)/plab(3) c vertice(2)=vertice(2)-dz*plab(2)/plab(3) c*** Assign values to some ntuple variables x0(1)=vertice(1) x0(2)=vertice(2) xp(1)=plab(1) xp(2)=plab(2) xp(3)=plab(3) if (sqrt(xp(1)**2+xp(2)**2+xp(3)**2).eq.0.) then write(6,*)' ***** Severe ERROR: momentum = 0.!! ',xp(1),xp(2),xp(3) write(6,*) plabmo, db, evathe, evaphi xp(3)=plabmo endif thetapar=acos(xp(3)/sqrt(xp(1)**2+xp(2)**2+xp(3)**2)) phipar=uatan2(xp(2),xp(1)) xc=XFREON(3)*tan(thetapar)*cos(phipar) yc=XFREON(3)*tan(thetapar)*sin(phipar) xpar=w*tan(thetapar)*cos(phipar) ypar=w*tan(thetapar)*sin(phipar) * write(6,*) '********* ',evathe,thetapar,evaphi,phipar,' ********' endif 1335 continue c*** in clas12 radiator acceptance ?? (TB improved) x_sec = vertice(1)-XOFREO y_sec = vertice(2) r_sec = SQRT(x_sec*x_sec+y_sec*y_sec) a_sec = ATAN2(y_sec,x_sec)*180./3.1415 * write(6,*) n_out, x_sec, y_sec, r_sec, a_sec ** try again if out of the clas acceptance if(ifine.eq.0) then IF ((ABS(a_sec) .GT. ANFREO) .OR. (r_sec .LT. R1FREO) .OR. (r_sec .GT. R2FREO)) THEN n_out = n_out + 1 if (n_out .LT. 100) then goto 3321 endif write (6,*) 'SEVERE WARNING: ',n_out, ' GUKINE EXTRACTIONS OUT OF ACCEPTANCE' ENDIF ENDIF c*** vertice iniziale e cinematica call gsvert (vertice,0,0,ub,1,nvtx) call gskine (plab,typpa,nvtx,ub,1,nt) if (iverbo.gt.1) then call gpvert(0) call gpkine(0) endif call gmedia(vertice,iiiee) ipar=typpa typpa=typpa+3 if (iverbo.gt.3) then write(6,*)' IVERBO-3 : gukine end' endif return 1136 continue ifine=1 WRITE(6,*)' INFO : End of input track File, program will end!' WRITE(6,*)' INFO : Last event analyzed:',ivent WRITE(6,*)' WARNING: LAST track must be skipped in the analysis!' c 111 format(6F10.4) goto 1335 end c************************************************************************** subroutine gutrev c************************************************************************** #undef CERNLIB_GEANT321_GCNUM_INC C #include "geant321/gcnum.inc" C NPART = 100 c write(6,*) 'SUBROUTINE GUTREV' call gtreve return end c*************************************************************************** subroutine gutrak c*************************************************************************** call gtracka c write(6,*) 'SUBROUTINE GUTRAK' return end c*************************************************************************** subroutine guskip (iskip) c*************************************************************************** iskip=0 c write(6,*) 'SUBROUTINE GUSKIP' return end c**************************************************************************** subroutine gustep c**************************************************************************** #undef CERNLIB_GEANT321_GCKING_INC #include "geant321/gcking.inc" #undef CERNLIB_GEANT321_GCTRAK_INC #include "geant321/gctrak.inc" #undef CERNLIB_GEANT321_GCKINE_INC #include "geant321/gckine.inc" #undef CERNLIB_GEANT321_GCTMED_INC #include "geant321/gctmed.inc" #undef CERNLIB_GEANT321_GCSETS_INC #include "geant321/gcsets.inc" integer ievafl integer typpa real xpp(3) real effqe real XPROXI(3),XFREON(3),XQUARZO(3),XGAP(3),XDET(3),XDRIFT(3),ROTRIC(6),XCERE(3) real R1FREO,R2FREO,XOFREO,ANFREO,R1DETE,R2DETE,XODETE,A1DETE,A2DETE common/matsiz/XPROXI,XFREON,XQUARZO,XGAP,XDET,XDRIFT,ROTRIC,XCERE common/clas12/R1FREO,R2FREO,XOFREO,ANFREO,R1DETE,R2DETE,XODETE,A1DETE,A2DETE common/inpdat/plabmi,dz,sigfax,sigfay,sigfte,sigfph,plabma,typpa,plabmo common/riccar/nrpho,ntpho,nnttph,w,iverbo common/qetab/eqe(14),qe(14) parameter (maxgamma=2000) common/trkpar/ivent,ipar,xp(3),x0(2),ngtot,nglost,ngamma,xg(maxgamma),yg(maxgamma), + thetag(maxgamma),egamma(maxgamma),z0g(maxgamma), + theta0g(maxgamma),nfreon,xpar,ypar,thetapar,phipar,xc,yc, + thetac(maxgamma), + phic(maxgamma),x(maxgamma),y(maxgamma),b(maxgamma), + etha(maxgamma) dimension xpho(3),vec(3),eff(maxgamma) logical pri pri=.false. if(ivent.eq.6)pri=.false. if(pri)write(*,*) 'SUBROUTINE GUSTEP',ivent,ngphot,ipart ievafl = 0 if (iverbo.gt.4) then write(6,*)' IVERBO-4 : gustep started' endif irifra=0 c*** elimina fotoni riflessi if (nmec .gt. 0) then do i=1,nmec if ((lmec(i).eq.102).or.(lmec(i).eq.106)) irifra=irifra+1 enddo if (irifra .gt. 1) then nrpho=nrpho+1 istop=1 ! blocca fotoni riflessi endif endif c*** store tracking steps call gsxyz c*** salva nello stack tutte le particelle generate durante il passo c*** ma non i fotoni cherenkov !! if (ngkine.gt.0) then if (iverbo.gt.1) write(*,*) ' ---> SECONDARIES!!! ' if (ngkine.lt.5) then call gsking(0) else istop=1 endif end if c*** salva nello stack i fotoni Cherenkov generati nel passo c*** (analogo a gsking ma per i Cherenkov !!) do j=1,3 xpp(j)=xp(j)/sqrt(xp(1)**2+xp(2)**2+xp(3)**2) enddo xedu=2.*XFREON(3) if ((ngphot.gt.0)) then ntpho=ntpho+ngphot do i=1,ngphot c*** photon moving to positive z-direction, and generated in the freon if ((xphot(6,i).gt.0.).and.(ngamma.lt.maxgamma).and.(xphot(3,i).le.xedu)) then ngamma=ngamma+1 ITRA=ngamma call gskpho(i) c*** assign values to some PAW variables (vertex) xpho(1)=xphot(4,i)/sqrt(xphot(4,i)**2+xphot(5,i)**2+ xphot(6,i)**2) xpho(2)=xphot(5,i)/sqrt(xphot(4,i)**2+xphot(5,i)**2+ xphot(6,i)**2) xpho(3)=xphot(6,i)/sqrt(xphot(4,i)**2+xphot(5,i)**2+ xphot(6,i)**2) z0g(ngamma)=xphot(3,i) theta0g(ngamma)=asin(sqrt((xpho(2)*xpp(1)-xpho(1)*xpp(2))**2+ + (xpho(1)*xpp(3)-xpho(3)*xpp(1))**2+ + (xpho(3)*xpp(2)-xpho(2)*xpp(3))**2))*180./3.1415 egamma(ngamma)=0. if(pri)write(*,1011)i,xphot(1,i),xphot(2,i),xphot(3,i),xphot(7,i)*1.e9,xpho(1),xpho(2),xpho(3),theta0g(ngamma) else nrpho=nrpho+1 if ((ngamma.ge.maxgamma) .and. (ievafl .eq. 0)) then write(6,*)' WARNING :',ivent/3,' maxgamma photon reached',ngphot, ' ', ngamma,' not stored' ievafl=1 endif endif enddo endif 1011 FORMAT('emesso ',i6,3f10.2,4f10.3,4f10.3) 1012 FORMAT('phodet ',3f10.2,2f10.3) 1013 FORMAT('trovato',3f10.2,2f10.3) c*** photon hits the detector plane if (ipart .eq. 50) then C write (6,*) 'Photon hit: ',ipart, isvol, iset, idet C write (6,*) vect(1), vect(2), vect(3) CEV if ((isvol .ne. 0) .and. ( iset .ne. 0 ) .and. (idet .ne. 0)) then xedu=2.*(XFREON(3)+XQUARZO(3)+XGAP(3))-0.0002 if (vect(3) .GT. xedu) THEN C is photon in clas12 detector phase space ?? C x_sec = vect(1) + XODETE C y_sec = vect(2) C r_sec = SQRT(x_sec*x_sec+y_sec*y_sec) C a_sec = ATAN2(y_sec,x_sec)*180./3.1415 C write (6,*) ' Test ',x_sec, y_sec, r_sec, a_sec C IF (a_sec .GE. A1DETE) .AND. (a_sec .LE. A2DETE) .AND. (r_sec .GE. R1DETE) .AND. (r_sec .LE. R2DETE)) THEN call ranlux(eff,1) if ((vect(7).gt.eqe(1)).and.(vect(7).lt.eqe(14))) then do j=2,14 if (vect(7).lt.eqe(j)) then effqe=qe(j-1)+(qe(j)-qe(j-1))/(eqe(j)-eqe(j-1))*(vect(7)-eqe(j-1)) go to 333 end if enddo endif C ENDIF effqe=0. 333 continue if(pri)write(*,1012)vect(1),vect(2),vect(3),vect(7)*1.e9,effqe if (eff(1).lt.effqe) then egamma(itra)=vect(7) ! energia del fotone xg(itra)=vect(1) ! hit coordinates yg(itra)=vect(2) if (sqrt(vect(4)**2+vect(5)**2+vect(6)**2).eq.0.) then write(6,*)' **** SEVERE ERROR, gamma momentum = 0. !! set to 1.' vect(6)=1.E-9 endif vec(1)=vect(4)/sqrt(vect(4)**2+vect(5)**2+vect(6)**2) vec(2)=vect(5)/sqrt(vect(4)**2+vect(5)**2+vect(6)**2) vec(3)=vect(6)/sqrt(vect(4)**2+vect(5)**2+vect(6)**2) thetag(itra)=asin(sqrt(vec(1)**2+vec(2)**2))*180./3.1415 if(pri)write(*,1013)vect(1),vect(2),vect(3),vect(7)*1.e9,effqe else egamma(itra)=-1. endif endif endif if (iverbo.gt.4) then write(6,*)' IVERBO-4 : gustep end' endif return end c************************************************************************** subroutine gudigi c************************************************************************** c write(6,*) 'SUBROUTINE GUDIGI' return end c************************************************************************** subroutine guout c************************************************************************** #undef CERNLIB_GEANT321_GCONST_INC C #include "geant321/gconst.inc" #undef CERNLIB_GEANT321_GCKINE_INC C #include "geant321/gckine.inc" #undef CERNLIB_GEANT321_GCFLAG_INC #include "geant321/gcflag.inc" character*25 str,opzione SAVE opzione integer typpa integer nngamma,iric integer ifine real*8 fricd1, fricd2 integer fricd3, fricd4 real*8 bric c*** indici di rifrazione medi real*8 xxnf, xxnf2, xxnq2, xxng2, xthmax real XPROXI(3),XFREON(3),XQUARZO(3),XGAP(3),XDET(3),XDRIFT(3),ROTRIC(6),XCERE(3) common/matsiz/XPROXI,XFREON,XQUARZO,XGAP,XDET,XDRIFT,ROTRIC,XCERE common /ric/bric common /xxnn/xxnf, xxnf2, xxnq2, xxng2, xthmax common/ancora/ifine common/inpdat/plabmi,dz,sigfax,sigfay,sigfte,sigfph,plabma,typpa,plabmo common/grafic/inter,kwtip common/view/teta,fi,psi,u0,v0,su,sv common/color/ICOLG common/riccar/nrpho,ntpho,nnttph,w,iverbo parameter (maxgamma=2000) common/trkpar/ivent,ipar,xp(3),x0(2),ngtot,nglost,ngamma,xg(maxgamma),yg(maxgamma), + thetag(maxgamma),egamma(maxgamma),z0g(maxgamma), + theta0g(maxgamma),nfreon,xpar,ypar,thetapar,phipar,xc,yc, + thetac(maxgamma), + phic(maxgamma),x(maxgamma),y(maxgamma),b(maxgamma), + etha(maxgamma) parameter (maxpad=1000) common/padpar/npad,ipad(maxpad),jpad(maxpad),xgp(maxpad),ygp(maxpad),ngpad(maxpad),rthpad(maxpad) real pad_nx, pad_ny real pad_dx, pad_dy,pad_gx,pad_gy,pad_sx,pad_sy common/padgeo/pad_dx,pad_dy,pad_gx,pad_gy,pad_sx,pad_sy,pad_nx,pad_ny external fric external ipadij save ntrial save iric c write(6,*) 'SUBROUTINE GUOUT' * first call, inter > 1 means number of events in background if ((iric .eq. 0).and.(inter .gt. 1)) iric=inter if (iverbo.gt.3) then write(6,*)' IVERBO-3 : guout started' endif call gfhsta c write(*,*)'##EVE##',ivent c*** Filling ntuples nngamma=0 nglost=0 nfreon=0 do i=1,ngamma C write (6,*) egamma(i), z0g(i), 2.*XFREON(3) if (egamma(i).gt.0) then if (z0g(i).le.(2.*XFREON(3))) then nfreon=nfreon+1 endif nngamma=nngamma+1 egamma(nngamma)=egamma(i) z0g(nngamma)=z0g(i) theta0g(nngamma)=theta0g(i) thetag(nngamma)=thetag(i) xg(nngamma)=xg(i) yg(nngamma)=yg(i) else if (egamma(i).eq.0) then nglost=nglost+1 endif endif enddo ngtot=ntpho ngamma=nngamma c*** digitizzazione detector npad=0 * write(*,*) 'NGAM: ',ngamma do i=1,ngamma ievai=ipadij(xg(i),yg(i),ig,jg,pad_dx,pad_dy,pad_gx,pad_gy,pad_sx,pad_sy, pad_nx, pad_ny) * write(*,*) ievai,xg(i),yg(i),pad_dx,pad_dy,pad_gx,pad_gy,pad_sx,pad_sy,ig,jg if (ievai.eq.1) then do j=1,npad if((ig.eq.ipad(j)).and.(jg.eq.jpad(j))) then ngpad(j)=ngpad(j)+1 go to 123 endif enddo npad=npad+1 ipad(npad)=ig jpad(npad)=jg ngpad(npad)=1 endif 123 continue enddo * write(*,*) 'NPAD: ',npad c******** angle reconstruction (all pe) do i=1,ngamma x(i)=xg(i)-x0(1) y(i)=yg(i)-x0(2) b(i)=sqrt((x(i)-xc)**2+(y(i)-yc)**2) phic(i)=ATAN2((y(i)-yc),(x(i)-xc)) bric=dble(b(i)) fricd1=0. fricd2=0.001 fricd3=10000 fricd4=1 thetac(i)=REAL(dzerox(fricd1,xthmax,fricd2,fricd3,fric,fricd4)) xdummy1=sin(thetapar)*sin(thetac(i))*cos(phic(i)-phipar) xdummy2=cos(thetapar)*cos(thetac(i)) etha(i)=acos(xdummy1+xdummy2) enddo c******** angle reconstruction (all PAD) do i=1,npad call padxy(ipad(i),jpad(i),xeva,yeva,pad_dx,pad_dy,pad_gx,pad_gy,pad_sx,pad_sy, pad_nx, pad_ny) * write(*,*) ipad(i),jpad(i),pad_dx,pad_dy,pad_gx,pad_gy,pad_sx,pad_sy,xeva,yeva xevax=xeva-x0(1) yevay=yeva-x0(2) xgp(i)=xeva ygp(i)=yeva bevab=sqrt((xevax-xc)**2+(yevay-yc)**2) pevap=ATAN2((yevay-yc),(xevax-xc)) bric=dble(bevab) fricd1=0. fricd2=0.001 fricd3=10000 fricd4=1 tevat=REAL(dzerox(fricd1,xthmax,fricd2,fricd3,fric,fricd4)) xdummy1=sin(thetapar)*sin(tevat)*cos(pevap-phipar) xdummy2=cos(thetapar)*cos(tevat) rthpad(i)=acos(xdummy1+xdummy2) enddo ntrial=ntrial+1 ivent=ntrial if(ifine.eq.0) call hfnt(10) ngamma=0 if (ifine.eq.1) then ieorun=1 else if (iric .gt. 0) iric=iric-1 if ((iric.eq.0).and.(inter .gt. 1)) ieorun=1 * if ((inter.eq.1).and.(ieorun.eq.0).and.(iric.eq.0)) then if ((ieorun.eq.0).and.(iric.eq.0)) then write(6,*) ' *** TRIAL ',ntrial, ' ********************' write(*,*) ' --- Total Ch photons =',ntpho, ' reflected =', nrpho c..... define viewport /window dimensione (see HIGZ manual) if (ntrial.eq.1) then c call ISWKVP(1,0.,1.,0.,1.) c call iswkwn(1,0.,1.,0.,1.) opzione='A' call ISWN(10,0.,20.,0.,20.) call ISVP(10,.005,.495,.005,.495) call ISWN(11,0.,20.,0.,20.) call ISVP(11,.505,.99,.005,.495) call ISWN(12,0.,20.,0.,20.) call ISVP(12,.005,.495,.505,.995) call ISWN(13,0.,20.,0.,20.) call ISVP(13,.505,.995,.505,.995) call ISWN(15,0.,20.,0.,20.) call ISVP(15,.005,.995,.005,.995) endif call gdopt ('TRAK','LINE') CALL GDOPT ('HIDE','ON') call gsatt ('PROX','SEEN',0) call gsatt ('FREO','COLO',ICOLG) call gsatt ('QUAR','COLO',ICOLG+1) call gsatt ('GAPP','SEEN',0) call gsatt ('DETE','COLO',ICOLG+3) CEV call gsatt ('WIRX','SEEN',0) CEV call gsatt ('WIRY','SEEN',0) CALL ISELNT(15) if ((opzione .eq. '1').or.(opzione.eq.'A')) then if (opzione .eq. 'A') then CALL ISELNT(13) endif call gdraw ('PROX',teta,fi,psi,u0,v0,su,sv) call gdxyz (0) call IGTEXT(19.,19.,'1',.5,0.,'L') endif if ((opzione .eq. '2').or.(opzione.eq.'A')) then if (opzione .eq. 'A') then CALL ISELNT(12) endif call gdraw ('PROX',90.,90.,270.,u0,v0,su,sv) call gdxyz (0) call gdscal (10.,1.) call IGTEXT(19.,19.,'2',.5,0.,'L') endif if ((opzione .eq. '3').or.(opzione.eq.'A')) then if (opzione .eq. 'A') then CALL ISELNT(11) endif call gdopt('TRAK','POIN') call gdraw ('PROX',teta,fi,psi,u0,v0,su,sv) call gdxyz (0) call IGTEXT(19.,19.,'3',.5,0.,'L') endif if ((opzione .eq. '4').or.(opzione.eq.'A')) then if (opzione .eq. 'A') then CALL ISELNT(10) endif call gdopt('TRAK','POIN') call gdraw ('PROX',0.,0.,0.,u0,v0,su,sv) call gdxyz (0) call IGTEXT(19.,19.,'4',.5,0.,'L') endif if (opzione .eq. 'Q') then CALL GDSPEC('PROX') write(*,*)' Last Run then Exit' ieorun = 1 endif if (opzione .eq. 'V') then write(*,*)' Change Verbose Level: (current is ',iverbo,') ' read(*,*)iverbo opzione='A' endif if (opzione .eq. 'B') then write(*,*)' Go in Background, give number of trigger then ' read(*,*)iric opzione='A' endif CALL ISELNT(15) call IGTEXT(0.1,2.2,'Menu"J# (+ RETURN)',.3,0.,'L') call IGTEXT(0.5,1.6, '1,2,3,4 - single view ',.3,0.,'L') call IGTEXT(0.5,1.1,'A - all, V - verbose level ',.3,0.,'L') call IGTEXT(0.5,.6, ' - current choice, B - background ',.3,0.,'L') call IGTEXT(0.5,0.1,'Q - last run then quit',.3,0.,'L') str=' ' call irqst (1,1,istat,nch,str) if (nch .gt. 0) then opzione=str endif write(6,*)' Opzione = ',opzione call iclrwk (1,0) end if end if if (iverbo.gt.3) then write(6,*)' IVERBO-3 : guout end' endif return end c************************************************************************* C Descrizione geometrica e proprieta' fisiche del rivelatore c************************************************************************* subroutine proxifoc #undef CERNLIB_GEANT321_GCTRAK_INC C #include "geant321/gctrak.inc" #undef CERNLIB_GEANT321_GCBANK_INC #include "geant321/gcbank.inc" integer typpa character*4 CHNMSV(1),CHNAMH(1) integer NBITSV(1) integer NBITSH(1) real ORIG(1),FACT(1) dimension ub(1) real afreon(2),zfreon(2),dens1(1),wfreon(2) real aquarzo(2),zquarzo(2),dens2(1),wquarzo(2) real ameta(2),zmeta(2),densmeta(1),wmeta(2) real XPROXI(3),XFREON(3),XQUARZO(3),XGAP(3),XDET(3),XDRIFT(3),XCERE(3) real R1FREO,R2FREO,XOFREO,ANFREO,R1DETE,R2DETE,XODETE,A1DETE,A2DETE real ROTRIC(6) c*** dim. parametri ottici dei volumi usati real ppfreo(14),abfreo(14),effreo(14),rifreo(14) real abquar(14),efquar(14),riquar(14) real abgap(14),efgap(14),rigap(14) real abdet(14),efdet(14),ridet(14) real xf_min, xf_max, dxf, dyf, xd_min,xd_max, dxd, dyd c*** indici di rifrazione medi real*8 xxnf, xxnf2, xxnq2, xxng2, xthmax real xdummy real psizex,psizey,pmaxx,pmaxy integer ipd,jpd,iii real mxp,myp real xwire(3), ywire(3), rbord(3) real pad_nx, pad_ny real pad_dx, pad_dy,pad_gx,pad_gy,pad_sx,pad_sy integer rad_ty common/radtyp/rad_ty common/inpdat/plabmi,dz,sigfax,sigfay,sigfte,sigfph,plabma,typpa,plabmo common/matsiz/XPROXI,XFREON,XQUARZO,XGAP,XDET,XDRIFT,ROTRIC,XCERE common/clas12/R1FREO,R2FREO,XOFREO,ANFREO,R1DETE,R2DETE,XODETE,A1DETE,A2DETE common /xxnn/xxnf, xxnf2, xxnq2, xxng2, xthmax common/padgeo/pad_dx,pad_dy,pad_gx,pad_gy,pad_sx,pad_sy,pad_nx,pad_ny c*** Inserimento dati data afreon,zfreon/12.011,18.998,6.000,9.000/ data dens1,wfreon/1.68,6.0,14.0/ data aquarzo,zquarzo/28.0855,15.9994,14.000,8.000/ data dens2,wquarzo/4.40,1.0,2.0/ data ameta,zmeta/12.011,1.00794,6.000,1.000/ data densmeta,wmeta/0.4241,1.0,4.0/ data CHNMSV/'DETE'/ data NBITSV/6/ data CHNAMH/'X '/ data NBITSH/16/ data ORIG/0./ data FACT/1./ write(6,*) 'SUBROUTINE PROXIFOC' c*** lettura parametri ottici dei materiali c Attenzione, l'indice di rifrazione del FREON e' ricalcolato sotto C ppfreo=photon energy in GeV in freon C abfreo=absorption lenght for freon (cm?... to be verified!) C effreo=efficiency for freon (always 1) C rifreo=refractive index for freon C abquar=absorption lenght for quartz (cm?... to be verified!) C efquar=efficiency for quartz (always 1) C riquar=refractive index for quartz C abgap=absorption lenght for gas (cm?... to be verified!) C efgap=efficiency for gas (always 1) C rigap=refractive index for gas write (*,*) 'Read optical Parameters from: optic.dat' open (unit=1,file='input/optic.dat',status='old') read(1,*) (ppfreo(i),abfreo(i),effreo(i),rifreo(i),abquar(i),efquar(i), + riquar(i),abgap(i),efgap(i),rigap(i),abdet(i), + efdet(i),ridet(i), i=1,14) close (unit=1) xedu=0. ! dummy value c*** valori medi indici di rifrazione: xxnf = 0. xxnq2 = 0. xxng2 = 0. do i=1,14 C formula di AntonelloDi Mauro / Settembre 2004 C Freon C6F14: if (rad_ty .EQ. 0) THEN C Fit dati Evaristo c rifreo(i) = 1.177 + 0.0172 * ppfreo(i) * 1.e9 C Fit dati Moyssdes rifreo(i)=1.2511+2.1565/(1240/(ppfreo(i) * 1.e9)-129.92) WRITE(*,*) 'Freon Type = C6F14 (Used in Hypernuclear Exp.)' ENDIF C Freon C5F12: IF (rad_ty .EQ. 1) THEN xdummy = 68.247 / ( 21.364**2 - (ppfreo(i) * 1.e9)**2) rifreo(i) = sqrt((1. + 2. * xdummy)/(1. - xdummy)) WRITE(*,*) 'Freon Type = C5F12 (to be cooled down)' ENDIF xxnf = xxnf + dble(rifreo(i)) xxnq2 = xxnq2 + dble(riquar(i)) xxng2 = xxng2 + dble(rigap(i)) enddo xxnf = xxnf / 14. xxnf2 = xxnf * xxnf xxnq2 = xxnq2 / 14. xxnq2 = xxnq2 * xxnq2 xxng2 = xxng2 / 14. xthmax = DASIN(xxng2/xxnf)-0.001 xxng2 = xxng2 * xxng2 write(6,*) ' ---- Average Refraction Indeces:' write(6,*) ' Freon = ',xxnf,', Quarzo = ',dsqrt(xxnq2), ', Gap = ',dsqrt(xxng2) write(6,*) ' Max Angle = ',xthmax,' rad' c*** definisce i materiali non standard c*** definizione del freon come mistura call gsmixt (100,'FREON$',afreon,zfreon,dens1,-2,wfreon) c*** definizione del quarzo come mistura call gsmixt (101,'QUARZO$',aquarzo,zquarzo,dens2,-2,wquarzo) c*** definizione del metano come mistura call gsmixt (102,'METANO$',ameta,zmeta,densmeta,-2,wmeta) c*** tracking medium parameters igauto=1 dmb=0.01 deb=0.7 epb=0.01 stb=0.01 ub(1)=0. call gstmed (1,'ARIA$', 15, 0,0,0.,0.,dmb,deb,epb,stb,ub,1) call gstmed (2,'VUOTO$', 16, 0,0,0.,0.,dmb,deb,epb,stb,ub,1) call gstmed (3,'FREON$', 100,0,0,0.,0.,dmb,deb,epb,stb,ub,1) call gstmed (4,'QUARZO$',101,0,0,0.,0.,dmb,deb,epb,stb,ub,1) call gstmed (5,'ALLU$', 9, 1,0,0.,0.,dmb,deb,epb,stb,ub,1) call gstmed (6,'META$', 102,0,0,0.,0.,dmb,deb,epb,stb,ub,1) c*** associa materiali e volumi c*** attenzione!!! Specificare il volume madre per primo! c*** ricalcola volume madre e resto sulla base dei parametri clas12 xf_min = XOFREO xf_max = XOFREO + R2FREO IF (ABS(xf_min) .GT. ABS(xf_max)) THEN dxf = ABS(xf_min) ELSE dxf = ABS(xf_max) ENDIF dyf = R2FREO * SIN(ANFREO/180.*3.1415) write (6,*) 'Freon: ',xf_min, xf_max, dxf, dyf xd_min = XODETE xd_max = XODETE + R2DETE IF (ABS(xd_min) .GT. ABS(xd_max)) THEN dxd = ABS(xd_min) ELSE dxd = ABS(xd_max) ENDIF dyd = R2DETE * (SIN(A2DETE/180.*3.1415)-SIN(A1DETE/180.*3.1415)) write (6,*) 'Detec: ',xd_min, xd_max, dxd, dyd XPROXI(3)=XFREON(3)+XQUARZO(3)+XGAP(3)+XDET(3)+10. XPROXI(1) = dxd+10. XPROXI(2) = dyd+10. XFREON(1) = dxf+1. XFREON(2) = dyf+1. XQUARZO(1) = dxf+2. XQUARZO(2) = dyf+2. XGAP(1) = dxd+1. XGAP(2) = dyd+1. pad_sx = dxd+0.1 pad_sy = dyd+0.1 pad_nx = 1 pad_ny = 1 rbord(1)=pad_sx rbord(2)=pad_sy rbord(3)=XDET(3) c*** wires xwire(1)=0.0 xwire(2)=0.01 xwire(3)=XGAP(1) xstep=.2 xdetdi=.4 ywire(1)=0.0 ! internal diameter ywire(2)=0.002 ! external diameter ywire(3)=XGAP(1) ! wire length ystep=.4 ! wire step ydetdi=.2 ! distance along track from the detector call gsvolu ('PROX','BOX ',2,XPROXI, 3,iproxi) !volume madre call gsvolu ('FREO','BOX ',3,XFREON, 3,ifreon) !radiatore di freon call gsvolu ('QUAR','BOX ',4,XQUARZO,3,iquarzo) !finestra in quarzo call gsvolu ('GAPP','BOX ',6,XGAP, 3,igap) !gap (metano) call gsvolu ('DETE','BOX ',5,rbord, 3,idete) !detector CEV call gsvolu ('WIRX', 'TUBE', 5, xwire, 3, iwirex) ! single wire cathode CEV call gsvolu ('WIRY', 'TUBE', 5, ywire, 3, iwirex) ! single wire anode C call gsvolu ('DRIF','BOX ',2,XDRIFT, 3,idrift) !drift + scintill. C call gsvolu ('CERE','BOX ',2,XCERE, 3,icere) !Gas Cerenkov c*** posiziona i volumi c write (*,*) 'Valori: ',ROTRIC(1),ROTRIC(2),ROTRIC(3),ROTRIC(4),ROTRIC(5),ROTRIC(6) call gsrotm (88,ROTRIC(1),ROTRIC(2),ROTRIC(3),ROTRIC(4),ROTRIC(5),ROTRIC(6)) call gsrotm (89, 180., 0., 90., 90., 90., 0.) call gsrotm (90, 90., 0., 180., 0., 90., 90.) xedu=XFREON(3) write(6,*)' z_freon center = ',xedu call gspos ('FREO',1,'PROX',0.,0.,xedu,0,'ONLY') xedu=2.*XFREON(3)+XQUARZO(3)-0.0001 write(6,*)' z_quartz center = ',xedu call gspos ('QUAR',1,'PROX',0.,0.,xedu,0,'ONLY') xedu=2.*(XFREON(3)+XQUARZO(3))+XGAP(3)-0.0001 write(6,*)' z_gap center = ',xedu call gspos ('GAPP',1,'PROX',0.,0.,xedu,0,'ONLY') C zedu= 2.*XFREON(3)+XCERE(3) .... frame laterali che fanno ombra sulle pad C xedu= XFREON(1)+56.+XCERE(1) C yedu= XFREON(2)+56.+XCERE(2) xedu=2.*(XFREON(3)+XQUARZO(3)+XGAP(3))+XDET(3)-0.0002 write(6,*)' z_detector (pad) center = ',xedu C POSITIONING OF PAD psizex=2.*pad_sx+pad_gx psizey=2.*pad_sy+pad_gy pmaxx=pad_nx*psizex/2. pmaxy=pad_ny*psizey/2. iii=1 do jpd=1,pad_ny myp = -pmaxy+(psizey)*(float(jpd)-0.5) do ipd=1,pad_nx mxp = -pmaxx+psizex*(float(ipd)-0.5) WRITE(*,*) 'Det ',iii,' xc: ',mxp,' yc: ',myp call gspos ('DETE',iii,'PROX',mxp,myp,xedu,0,'ONLY') iii=iii+1 enddo enddo C if (pad_gy .gt. 0) then C call gspos ('DETE',1,'PROX',0.,pad_sy+pad_gy/2.,xedu,0,'ONLY') C call gspos ('DETE',2,'PROX',2.*pad_sx+pad_gx,pad_sy+pad_gy/2.,xedu,0,'ONLY') C call gspos ('DETE',3,'PROX',-2.*pad_sx-pad_gx,pad_sy+pad_gy/2.,xedu,0,'ONLY') C call gspos ('DETE',4,'PROX',0.,-pad_sy-pad_gy/2.,xedu,0,'ONLY') C call gspos ('DETE',5,'PROX',2.*pad_sx+pad_gx,-pad_sy-pad_gy/2.,xedu,0,'ONLY') C call gspos ('DETE',6,'PROX',-2.*pad_sx-pad_gx,-pad_sy-pad_gy/2.,xedu,0,'ONLY') C else C if (pad_sy .gt. pad_sx) then C call gspos ('DETE',1,'PROX',-3.*pad_sx-1.5*pad_gx,0.,xedu,0,'ONLY') C call gspos ('DETE',2,'PROX',-pad_sx-0.5*pad_gx,0.,xedu,0,'ONLY') C call gspos ('DETE',3,'PROX',pad_sx+.5*pad_gx,0.,xedu,0,'ONLY') C call gspos ('DETE',4,'PROX',3.*pad_sx+1.5*pad_gx,0.,xedu,0,'ONLY') C else C call gspos ('DETE',1,'PROX',0.,0.,xedu,0,'ONLY') C call gspos ('DETE',2,'PROX',2.*pad_sx+pad_gx,0.,xedu,0,'ONLY') C call gspos ('DETE',3,'PROX',-2.*pad_sx-pad_gx,0.,xedu,0,'ONLY') C endif C endif CEV CEV write (*,*) 'Number of Wires: n_x = ',2.*XGAP(2)/xstep,', n_y = ',2.*XGAP(1)/ystep CEV CEV xedu = XGAP(3)-xdetdi CEV write (6,*)' x_Cathode Plane position: ',xedu CEV CEV do i=1,XGAP(2)*2./xstep CEV yoff = -XGAP(2) + (REAL(i) * xstep) - xstep/2. CEV call gspos ('WIRX', i, 'GAPP', 0., yoff, xedu, 89, 'ONLY') CEV enddo CEV xedu =XGAP(3)-ydetdi CEV write (6,*)' x_Anode Plane position: ',xedu CEV do i=1,XGAP(2)*2./ystep CEV xoff = -XGAP(2) + (REAL(i) * ystep) - ystep/2. CEV call gspos ('WIRY', i, 'GAPP', 0., xoff, xedu, 89, 'ONLY') CEV enddo C call gspos ('DRIF',1,'PROX',0.,10.,-90.,88,'ONLY') C call gspos ('DRIF',2,'PROX',0.,10.,-115.,88,'ONLY') C call gspos ('DRIF',3,'PROX',0.,0.,-6.,0,'ONLY') C call gspos ('DRIF',4,'PROX',0.,0.,73.,0,'ONLY') C call gspos ('CERE',1,'PROX',0.,0.,50.,0,'ONLY') c write (*,*) 'Matrice di rotazione: ' c call gprotm (88) c call irqst (1,1,istat,nch,str) c call iclrwk (0,kfl) c*** scrive in un file i parametri utili per pattern recognition c open(unit=32,name='proxi.par',type='new') c write(32,*) rw,qw,tgap,1.28,1.55,1.,deltaz c close(unit=32) c*** assegna materiale sensibile call gsdet('DESE','DETE',1,CHNMSV,NBITSV,1,100,100,iset,idet) call gsdeth('DESE','DETE',1,CHNAMH,NBITSH,ORIG,FACT) c*** chiusura descrizione geometrica call ggclos c*** dichiara proprieta' ottiche dei mezzi in uso C ab = absorbtion lenght C ef = efficiency C ri = refractive index call gsckov (3,14,ppfreo,abfreo,effreo,rifreo) call gsckov (4,14,ppfreo,abquar,efquar,riquar) call gsckov (6,14,ppfreo,abgap,efgap,rigap) call gsckov (5,14,ppfreo,abdet,efdet,ridet) c*** parametri in output call gpvolu (0) call gptmed (0) call gppart (0) call gpmate (100) call gpmate (101) call gpmate (102) call gpmate (9) return end c************************************************************************* c BYE c*************************************************************************