From Peter_Link at physik.tu-muenchen.de Tue Apr 6 11:09:06 1999 From: Peter_Link at physik.tu-muenchen.de (Peter Link) Date: Tue, 06 Apr 1999 11:09:06 +0200 Subject: V-selector Message-ID: <3709CF32.BB5EAA0C@physik.tu-muenchen.de> Dear Kristian, receiving the new release 1.1 of mcstas I was happy to find the new V-selector component. As you remember I sent you my own component definition for that in march,99. Looking at your solution I do not really see, why you use an analytical formula ( with approximations) for the transmission instead a MC-choice of the phase angle at the entrance position of the selector. To my feeling the advantage of the latter is that there's no further approximation neccessary ( apart from having a really black absorber on the blades). On the other hand I think there's no measureable change in the result, ?? As I'm looking to build more complicated versions of velocity selectors I would be glad to discuss this difference with you. Best regards Peter -- ****************************************** Dr. P. Link IPC Uni G?ttingen Au?enstelle Garching Neue Forschungs- Neutronenquelle Garching ZBE- FRM II- Bau 85747 Garching Tel. 089 2891 4622 Fax. 089 2895 4622 mailto: plink at physik.tu-muenchen.de ****************************************** From kristian.nielsen at risoe.dk Thu Apr 15 11:09:33 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: 15 Apr 1999 11:09:33 +0200 Subject: MCSTAS In-Reply-To: (message from Fabrizio Semadeni on Mon, 12 Apr 1999 09:02:12 +0100) Message-ID: > Since I'm not an expert in Perl, I prefer that you take a look from > yourself at > the possible incompatibilities, bugs and so on that occure here. You got a Ok, I looked at the problem. It seems that you may have a version problem with your perl installation. You have two versions installed: /progs/gnu/bin/perl (version 5.004) /progs/perl (version 5.005) The PGPLOT interface does not work with /progs/gnu/bin/perl, but it does work with /progs/perl. My guess is that PGPLOT was installed for the newer version and hence does not work with the old one. On my account, and hence presumably other accounts as well, the /progs/gnu/bin/perl version is first in the path, so when McStas is installed it chooses this version. I think that if /procs/perl is first in the path then a reinstall of McStas will fix the problems. I was able to succesfully run the graphics version of McStas using the correct perl version. Incidentally, the gcc compiler in /progs/gnu/bin also seems to have problems. Please ask again if you have any questions or the problem persists. - Kristian. From Ralph_Gilles at Physik.TU-Muenchen.DE Fri Apr 23 14:32:40 1999 From: Ralph_Gilles at Physik.TU-Muenchen.DE (Dr. Ralph Gilles) Date: Fri, 23 Apr 1999 14:32:40 +0200 Subject: supermirrors Message-ID: <37206868.108B8CA2@ph.tum.de> Dear Kristian, I have a question concerning the input of supermirrors. For example: R0 = 1.0 mh = 3, mv = 0 I made a few variations of these parameters to see the meaning of this example in detail. mv = 0 has still reflectivity because of R0 = 1.0. If I would like to calculate a supermirror with mh = 3 and a non reflecting material mv for the sidewalls (it make sense to check the contribution of the supermirror top and botton), I have the problem to set R0 = 0 for mv to really simulate this case. I know it is connected to the fall-off of reflectivity which is still by using m =0. The influence is of course not very big but in special cicumtances it could be very helpful. greetings, Ralph Gilles From kristian.nielsen at risoe.dk Fri Apr 23 15:59:43 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: 23 Apr 1999 15:59:43 +0200 Subject: supermirrors In-Reply-To: <37206868.108B8CA2@ph.tum.de> (Ralph_Gilles@Physik.TU-Muenchen.DE) Message-ID: Hi, > example in detail. mv = 0 has still reflectivity because of R0 = 1.0. If > I would like to calculate a supermirror with mh = 3 and a non reflecting > material mv for the sidewalls (it make sense to check the contribution For this, you should modify the Guide2 component to always absorb when m=0 (it is a very simple change). I have included some untested code below, could you try it and let me know if it works (so I can put it on the web page)? The only change is the added lines if((i <= 2 && mv == 0) || (i > 2 && mh == 0)) ABSORB; - Kristian. ------------------------------------------------------------------------ DEFINE COMPONENT Guide2 DEFINITION PARAMETERS (w1, h1, w2, h2, l, R0, Qc, alpha, mh, mv, W) SETTING PARAMETERS () STATE PARAMETERS (x,y,z,vx,vy,vz,t,s1,s2,p) TRACE %{ double t1,t2; /* Intersection times. */ double av,ah,bv,bh,cv1,cv2,ch1,ch2,d; /* Intermediate values */ double vdotn_v1,vdotn_v2,vdotn_h1,vdotn_h2; /* Dot products. */ int i; /* Which mirror hit? */ double q; /* Q [1/AA] of reflection */ double vlen2,nlen2; /* Vector lengths squared */ /* ToDo: These could be precalculated. */ double ww = .5*(w2 - w1), hh = .5*(h2 - h1); double whalf = .5*w1, hhalf = .5*h1; double lwhalf = l*whalf, lhhalf = l*hhalf; /* Propagate neutron to guide entrance. */ PROP_Z0; if(x <= -whalf || x >= whalf || y <= -hhalf || y >= hhalf) ABSORB; for(;;) { /* Compute the dot products of v and n for the four mirrors. */ av = l*vx; bv = ww*vz; ah = l*vy; bh = hh*vz; vdotn_v1 = bv + av; /* Left vertical */ vdotn_v2 = bv - av; /* Right vertical */ vdotn_h1 = bh + ah; /* Lower horizontal */ vdotn_h2 = bh - ah; /* Upper horizontal */ /* Compute the dot products of (O - r) and n as c1+c2 and c1-c2 */ cv1 = -whalf*l - z*ww; cv2 = x*l; ch1 = -hhalf*l - z*hh; ch2 = y*l; /* Compute intersection times. */ t1 = (l - z)/vz; i = 0; if(vdotn_v1 < 0 && (t2 = (cv1 - cv2)/vdotn_v1) < t1) { t1 = t2; i = 1; } if(vdotn_v2 < 0 && (t2 = (cv1 + cv2)/vdotn_v2) < t1) { t1 = t2; i = 2; } if(vdotn_h1 < 0 && (t2 = (ch1 - ch2)/vdotn_h1) < t1) { t1 = t2; i = 3; } if(vdotn_h2 < 0 && (t2 = (ch1 + ch2)/vdotn_h2) < t1) { t1 = t2; i = 4; } if(i == 0) break; /* Neutron left guide. */ PROP_DT(t1); switch(i) { case 1: /* Left vertical mirror */ nlen2 = l*l + ww*ww; q = V2Q*(-2)*vdotn_v1/sqrt(nlen2); d = 2*vdotn_v1/nlen2; vx = vx - d*l; vz = vz - d*ww; break; case 2: /* Right vertical mirror */ nlen2 = l*l + ww*ww; q = V2Q*(-2)*vdotn_v2/sqrt(nlen2); d = 2*vdotn_v2/nlen2; vx = vx + d*l; vz = vz - d*ww; break; case 3: /* Lower horizontal mirror */ nlen2 = l*l + hh*hh; q = V2Q*(-2)*vdotn_h1/sqrt(nlen2); d = 2*vdotn_h1/nlen2; vy = vy - d*l; vz = vz - d*hh; break; case 4: /* Upper horizontal mirror */ nlen2 = l*l + hh*hh; q = V2Q*(-2)*vdotn_h2/sqrt(nlen2); d = 2*vdotn_h2/nlen2; vy = vy + d*l; vz = vz - d*hh; break; } /* Now compute reflectivity. */ if((i <= 2 && mv == 0) || (i > 2 && mh == 0)) ABSORB; if(q > Qc) { double arg = (q - (i <= 2 ? mv : mh)*Qc)/W; if(arg < 10) p *= .5*(1-tanh(arg))*(1-alpha*(q-Qc)); else ABSORB; /* Cutoff ~ 1E-10 */ } p *= R0; } %} END From Georg_Artus at Physik.TU-Muenchen.DE Tue Apr 27 13:37:36 1999 From: Georg_Artus at Physik.TU-Muenchen.DE (Dr. Georg Artus) Date: Tue, 27 Apr 1999 11:37:36 +0000 Subject: mcplot Message-ID: <3725A180.8D3C8DBA@ph.tum.de> Hello Kristian, I have just installed version 1.1 of your program and I was delighted about the new mcplot routine! It works fine, but creating a *.ps file with 'P' ends up with: Click plot for full-screen view. Type 'P' (in graphics window) for hardcopy, 'Q' to quit. %PGPLOT, Unrecognized device type %PGPLOT, Invalid device specification: mcstas.ps/cps PGOPEN failed! at /usr/local/lib/mcstas/mcplotlib.pl line 120. Also the full screen display fails. What does 'Click plot ...' mean? Should I have a plot button somewhere? Do I have an incomplete installation of PGPLOT or PERL or what else might be the problem? Best wishes, Georg -- ********************************************* Dr. Georg Artus Technische Universitaet Muenchen Neue Forschungsneutronenquelle Garching D-85747 Garching Tel: +49 (0)89/289-14675 Fax: +49 (0)89/289-14666 or 12112 E-mail: gartus at ph.tum.de From kristian.nielsen at risoe.dk Tue Apr 27 14:02:03 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: 27 Apr 1999 14:02:03 +0200 Subject: mcplot In-Reply-To: <3725A180.8D3C8DBA@ph.tum.de> (Georg_Artus@Physik.TU-Muenchen.DE) Message-ID: > I have just installed version 1.1 of your program and I was delighted > about the new mcplot routine! It works fine, but creating a *.ps file > with 'P' ends up with: > > Click plot for full-screen view. > Type 'P' (in graphics window) for hardcopy, 'Q' to quit. > %PGPLOT, Unrecognized device type > %PGPLOT, Invalid device specification: mcstas.ps/cps > PGOPEN failed! at /usr/local/lib/mcstas/mcplotlib.pl line 120. This sounds very much like your pgplot was installed without color postscript support. Note this is the pgplot library itself, not the pgperl interface. Though this is somewhat strange, given that postscript is enabled by default in pgplot. Anyway, if you can find the file "drivers.list" somewhere in the pgplot directory, you can see if the lines mentioning color postscript is commented out: PSDRIV 3 /CPS PostScript printers, color, landscape Std F77 If the line starts with an exclamation sign ("!"), then color postscript is disabled, and pgplot must be re-installed with it enabled. (Erm ... wasn't it me who installed pgplot for you? I hope I did not make this mistake, that would be kind of embarrassing. Or maybe I only installed the pgperl interface, can't remember anymore.) > Also the full screen display fails. What does 'Click plot ...' mean? > Should I have a plot button somewhere? :-) "Click with the mouse button on the plot in the window ..." - I do not know where I got the telegram style from. - Kristian. From Georg_Artus at Physik.TU-Muenchen.DE Tue Apr 27 14:53:59 1999 From: Georg_Artus at Physik.TU-Muenchen.DE (Dr. Georg Artus) Date: Tue, 27 Apr 1999 12:53:59 +0000 Subject: mcplot References: <01JAJ3QDH2ZM921HZF@risoe.dk> Message-ID: <3725B367.3F4232DB@ph.tum.de> Hello Kristian, You are right, PSDRIV 3 /CPS was commented out. I'm wondering why you use color postscript as default output. I think most of the users only can use greyscale output having no color printer. So do we. Is it somehow possible to switch to greyscale output already on the screen? Otherwise I fear, that the dark red and the dark blue (opposite ends of the intensity scale) color will give very similar greyscales sent to a greyscale printer! > :-) "Click with the mouse button on the plot in the window ..." - I do > not know where I got the telegram style from. I tried that also, nothing happens with the window size. The latter is not very important, I was just playing. Georg -- ********************************************* Dr. Georg Artus Technische Universitaet Muenchen Neue Forschungsneutronenquelle Garching D-85747 Garching Tel: +49 (0)89/289-14675 Fax: +49 (0)89/289-14666 or 12112 E-mail: gartus at ph.tum.de From Georg_Artus at Physik.TU-Muenchen.DE Tue Apr 27 16:05:29 1999 From: Georg_Artus at Physik.TU-Muenchen.DE (Dr. Georg Artus) Date: Tue, 27 Apr 1999 14:05:29 +0000 Subject: mcplot / pgplot reinstalled References: <01JAJ3QDH2ZM921HZF@risoe.dk> Message-ID: <3725C429.BC592FC1@ph.tum.de> Hello Kristian, I have just reinstalled PGPLOT. Creating a ps file works now. But as I feared low and high intensity cannot be distinguished in the greyscale printout. Georg -- ********************************************* Dr. Georg Artus Technische Universitaet Muenchen Neue Forschungsneutronenquelle Garching D-85747 Garching Tel: +49 (0)89/289-14675 Fax: +49 (0)89/289-14666 or 12112 E-mail: gartus at ph.tum.de From Georg_Artus at Physik.TU-Muenchen.DE Wed Apr 28 09:51:39 1999 From: Georg_Artus at Physik.TU-Muenchen.DE (Dr. Georg Artus) Date: Wed, 28 Apr 1999 07:51:39 +0000 Subject: full screen problem ok References: <01JAJ3QDH2ZM921HZF@risoe.dk> Message-ID: <3726BE0B.52532F32@ph.tum.de> Hello Kristian, the full screen problem with mcplot is clear now! Yesterday I had only one detector in my instrument. So mcplot plotted this detector already in enlarged size! Therefor I thought clicking would enlarge the window to full screen size! I'm sorry for the confusion about this! It works if you have more than one detector in the instrument. Two more remarks: I had five detectors, so three in the upper row and two (left and middle position) in the lower. Clicking on the detecotr in the upper right corner (third one) produced a full window size of the detector in the middle of the lower row. Enlarging the third detector was not possible. I would like to know how one can use mcplot with a gscan output. I've tried several things but nothing works. I haven't found it in the manual. Best wishes, Georg -- ********************************************* Dr. Georg Artus Technische Universitaet Muenchen Neue Forschungsneutronenquelle Garching D-85747 Garching Tel: +49 (0)89/289-14675 Fax: +49 (0)89/289-14666 or 12112 E-mail: gartus at ph.tum.de From kristian.nielsen at risoe.dk Wed Apr 28 11:26:19 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: 28 Apr 1999 11:26:19 +0200 Subject: full screen problem ok In-Reply-To: <3726BE0B.52532F32@ph.tum.de> (Georg_Artus@Physik.TU-Muenchen.DE) Message-ID: > Two more remarks: > I had five detectors, so three in the upper row and two (left and middle > position) in the lower. Clicking on the detecotr in the upper right > corner (third one) produced a full window size of the detector in the > middle of the lower row. Enlarging the third detector was not possible. Ok, thanks for the bug report. I found a typo in the code that caused this, and fixed it. I also changed the default postscript output to B/W as you suggested, and improved the confusing message (all of them were simple to change). Thanks for your comments! If you want, you can now obtain updated files from the "known bugs" web page: http://elu-alf-2.risoe.dk/~elu-krni/mcstas/known-bugs.html > I would like to know how one can use mcplot with a gscan output. I've > tried several things but nothing works. I haven't found it in the > manual. Unfortunately, that is not currently possible. The gscan program is still a left-over from the very early McStas development, and it badly needs a major update, including an interface to mcplot. That will probably have to wait for the next major release of McStas, maybe later this year? - Kristian. From Daniel_Berger at physik.tu-muenchen.de Fri Apr 30 17:28:51 1999 From: Daniel_Berger at physik.tu-muenchen.de (Daniel Berger) Date: Fri, 30 Apr 1999 17:28:51 +0200 Subject: Windows Message-ID: <3729CC33.18297B15@physik.tu-muenchen.de> Dear Mr. Nielsen, I have got a PC with WINDOWS '95 and a Pentium-Processor. I have been trying to install McStas, but it wasn't sucessfull. Could you give me an exact description for installing MS under WINDOWS '95 and the software-/hardwarerequirements. Futhermore I don't understand the informations on your Homepage ("Installing on non-Unix systems). Thanks for helping D. Berger From Georg_Artus at Physik.TU-Muenchen.DE Thu Apr 1 11:22:03 1999 From: Georg_Artus at Physik.TU-Muenchen.DE (Dr. Georg Artus) Date: Thu, 01 Apr 1999 09:22:03 +0000 Subject: New release v1.1 of McStas References: <01J9HQK48M1Q90RMZK@risoe.dk> Message-ID: <37033ABB.792357DF@ph.tum.de> Hello Kristian, Thanks for the news about McStas! I have just finished the simulations with McStas! Our instrument is well defined now (as we think...). Again let me thank you for helping! McStas has been a very important feature in our work. The instrument and of course the simulation results will be presented at the N99 in Berlin and at the ECNS99 in Budapest. I'm sure you will be in Budapest. How about Berlin? Happy Easter Weekend! Georg -- ********************************************* Dr. Georg Artus Technische Universitaet Muenchen Neue Forschungsneutronenquelle Garching D-85747 Garching Tel: +49 (0)89/289-14675 Fax: +49 (0)89/289-14666 or 12112 E-mail: gartus at ph.tum.de From Ralph_Gilles at Physik.TU-Muenchen.DE Fri Apr 23 17:11:58 1999 From: Ralph_Gilles at Physik.TU-Muenchen.DE (Dr. Ralph Gilles) Date: Fri, 23 Apr 1999 17:11:58 +0200 Subject: supermirrors References: <01JADMOTBOB2920DA5@risoe.dk> Message-ID: <37208DBE.5329E8BA@ph.tum.de> Hello Kristian, thanks a lot. It is working very well. Georg and I tested the guide and it works very perfect. But in the input file first line: DEFINE COMPONENT GUIDE"3" instead of "2" is necessary. :-) greetings, Ralph > > example in detail. mv = 0 has still reflectivity because of R0 = 1.0. If > > I would like to calculate a supermirror with mh = 3 and a non reflecting > > material mv for the sidewalls (it make sense to check the contribution > > For this, you should modify the Guide2 component to always absorb when > m=0 (it is a very simple change). I have included some untested code > below, could you try it and let me know if it works (so I can put it on > the web page)? > > The only change is the added lines > > if((i <= 2 && mv == 0) || (i > 2 && mh == 0)) > ABSORB; > > - Kristian. > From kristian.nielsen at risoe.dk Mon May 3 08:21:56 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: 3 May 1999 08:21:56 +0200 Subject: Windows In-Reply-To: <3729CC33.18297B15@physik.tu-muenchen.de> (message from Daniel Berger on Fri, 30 Apr 1999 17:28:51 +0200) Message-ID: > I have got a PC with WINDOWS '95 and a Pentium-Processor. I have been > trying to install McStas, but it wasn't sucessfull. Could you give me > an exact description for installing MS under WINDOWS '95 and the > software-/hardwarerequirements. Futhermore I don't understand the > informations on your Homepage ("Installing on non-Unix systems). > Thanks for helping Unfortunately, you are mostly on your own here. We do not use Windows for McStas simulations at Ris? and cannot easily support that platform because of its poor compatibility. Users outside of Ris? have successfully compiled and run McStas on Windows (using Microsoft Visual C++ I believe). You can check the mailing list for details: http://neutron.risoe.dk/neutron-mc/ You should also try to ask for help on the mailing list (send mail to neutron-mc at risoe.dk). Alternatively, you could install Linux on the PC. For example, Debian GNU/Linux (http://www.debian.org/) works very well with McStas. The only requirement (for windows or any other platform) is a standard ANSI C compiler (GCC or Visual C++ should work fine). - Kristian. From kristian.nielsen at risoe.dk Mon May 3 10:06:58 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: 3 May 1999 10:06:58 +0200 Subject: Phase Space Transformation In-Reply-To: <37285439.2AC21BD@fz-juelich.de> (message from Oliver Kirstein on Thu, 29 Apr 1999 14:44:41 +0200) Message-ID: > Dear Kristian! > > It's been some time ago since Grenoble but I found some time and tried > to use the phase space transformation (which is nothing than bragg > scattering at a moving crystal with a large mosaic). It is based upon > Abd now the question: Can you crosscheck this routine if you have some > time to do it? I'm not sure if this would be interesting for anybody, Thanks for sending me your work on the phase space transformer component. I am sure that this would be useful for others when it is finished, especially if you find the time to document it. Ok, I looked at your code. I haven't completely understood the code, but I did manage to find a few problems, and by fixing them get some counts in the detector. The following points should get you going again. - When you compute mosaic_h, you use asin(Q/ki_betrag). This should be asin(Q/2.0/ki_betrag), should it not? - There is a problem in your code in that it makes assumptions about from which direction the neutron is coming (I think). This is partly inherited from the admittedly not particularly good Monochromator.comp code. For example, in your test instrument the neutron is actually intersecting the crystal from behind (v_x > 0), so when you compute G, you have to make it negative: gz_h = -Q * sin(mosaic_h); gx_h = -Q * cos(mosaic_h); - I do not understand your formula for computing kf. It seems to me that your code uses kf = 2*g - ki, that should be kf = 2*g + ki, I think(?). - In your instrument, the crystal is at the angle theta=69.2. If the crystal was stationary, the neutron would be scattered at an angle of 2*theta, with phase transformation it might be something else, but the detector should be at some other angle than theta, shouldn't it? A good way to test this is to put a 4PI PSD detector around the PST component. I put up my modifications to your code on the web server at http://neutron.risoe.dk/mcstas/support/kirstein/ but the code still needs more modifications. I hope this helps, please feel free to ask again if you have further questions. - Kristian. From kristian.nielsen at risoe.dk Fri May 7 10:36:12 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: 7 May 1999 10:36:12 +0200 Subject: velocity selector In-Reply-To: <36F74FD1.796AEE72@physik.tu-muenchen.de> (message from Peter Link on Tue, 23 Mar 1999 09:24:49 +0100) Message-ID: Dear Peter Link, Thank you again for your monochromator and velocity selector components. I now finally found the time to make them available on the McStas web page; hopefully they will be of use to others. The components were delayed because I was busy with the new release, but I will try to be quicker with submitted components in the future. - Kristian. From kadowaki at comp.metro-u.ac.jp Sun May 9 23:39:46 1999 From: kadowaki at comp.metro-u.ac.jp (kadowaki at comp.metro-u.ac.jp) Date: Mon, 10 May 1999 06:39:46 +0900 Subject: thank you Message-ID: <199905092139.GAA44586@comp.metro-u.ac.jp> Hello Kristian, Thank you for your solving the problem with our HP-UX computer. I am using McStas by IBM-AIX, PC-Linux, and DEC-Unix machines, which work completely in the same way and very well. Now I have another choice. We are using your remarkable McStas simulation soft to solve problems on Japanese version of ESS and to decide what is the best parameters of the to-be-upgraded thermal guide of our research reactor JRR3M. I like McStas. Sincerely Hiroaki Kadowaki From Daniel_Berger at Physik.TU-Muenchen.DE Mon Jun 14 18:52:06 1999 From: Daniel_Berger at Physik.TU-Muenchen.DE (Daniel Berger) Date: Mon, 14 Jun 1999 18:52:06 +0200 Subject: flux-sources Message-ID: <37653336.26A3DBB2@ph.tum.de> Hello, I'm working with the flux source at the moment and I have some questions: - What is the meaning of the parameter mccount in the Source_flux.comp - file? - I would like to calculate fluxes with the Source_div_file. How can I realize it? Thanks for answering these questions. Yours sincerely Daniel Berger From kristian.nielsen at risoe.dk Wed Jun 16 08:45:53 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: Wed, 16 Jun 1999 08:45:53 +0200 Subject: flux-sources In-Reply-To: <37653336.26A3DBB2@ph.tum.de> (message from Daniel Berger on Mon, 14 Jun 1999 18:52:06 +0200) Message-ID: <01JCGN9BRUPS94NEAQ@risoe.dk> > Date: Mon, 14 Jun 1999 18:52:06 +0200 > From: Daniel Berger > - What is the meaning of the parameter mccount in the Source_flux.comp - > file? This is an old (and undocumented) way to get the number of neutron histories to simulate (the value specified on the command line as -n or --ncount). From McStas version 1.1, this should be written as mcget_ncount() (as explained in the manual). > - I would like to calculate fluxes with the Source_div_file. How can I > realize it? (I assume you are referring to the "Source_div" component in the standard component library in McStas 1.1? Or are you using another component that I do not know about?) The component Source_div unfortunately does not correctly handle absolute flux. Is there any reason why you cannot use the Source_flux instead? One option is to put a Monitor_flux component just after the Source_div to measure the actual flux from the source, and then adjust the intensities appropriately. It would also be possible to come up with a formula for the flux from the Source_div component, and then modify the component to handle flux correctly just like the Source_flux component. If you want I will be happy to help you with this. - Kristian. From Daniel_Berger at Physik.TU-Muenchen.DE Wed Jun 16 16:52:36 1999 From: Daniel_Berger at Physik.TU-Muenchen.DE (Daniel Berger) Date: Wed, 16 Jun 1999 16:52:36 +0200 Subject: flux-sources References: <01JCGN9BRUPS94NEAQ@risoe.dk> Message-ID: <01JCJOOEOJNO94N8ZT@risoe.dk> Dear Mr. Nielsen, I have some problems to understand Source_flux. How does the programme calculate d_Phi/d_Lambda*d_Omega and how is d_Omega calculated? A further question is the influx of the source_area to d_Omega. Thanks for answering these questions. Yours sincerely Daniel Berger From kristian.nielsen at risoe.dk Fri Jun 18 13:24:53 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: Fri, 18 Jun 1999 13:24:53 +0200 Subject: flux-sources In-Reply-To: <3767BA33.8DF5FABF@ph.tum.de> (message from Daniel Berger on Wed, 16 Jun 1999 16:52:36 +0200) Message-ID: <01JCJPKXZQ7C94NTAD@risoe.dk> > Date: Wed, 16 Jun 1999 16:52:36 +0200 > From: Daniel Berger > I have some problems to understand Source_flux. > How does the programme calculate d_Phi/d_Lambda*d_Omega and how is d_Omega > calculated? A further question is the influx of the source_area to d_Omega. The basic idea in the Source_flux component is to start off mcget_ncount() neutron histories within a given wavelength range and solid angle. The initial neutron weight is then calculated so that the simulated intensity matches a given desired flux (neutrons per lambda per solid angle per source area). The solid angle (I think this is what you refer to with d_Omega?) is defined by a rectangle of size xw x yh at distance dist, so the maximum horizontal and vertical divergence may be computed as hdiv = atan(xw/(2.0*dist)) vdiv = atan(yh/(2.0*dist)) The solid angle within which neutron histories are emitted is therefore (4*hdiv*vdiv) (this implicitly makes approximations based on the assumptions that the angles are small and that the source dimensions are small compared to the distance dist. These approximation should probably be removed in a later version). The wavelength range delta_lambda of emitted neutron histories is easily computed from E0 and dE, as is the area of the source source_area from the source radius. Now the number of neutron histories per solid angle per wavelength per source area emitted in the simulation is F = mcget_ncount()/((4*hdiv*vdiv)*delta_lambda*source_area The simulated intensity in the simulation is p_in*F, where p_in is the initial weight of the neutron. Thus, to make the simulated intensity match a given flux, we should have p_in = flux/F This is what happens in the Source_flux component. I hope this answers your questions, otherwise feel free to ask again about details. - Kristian. From o.kirstein at fz-juelich.de Wed Jun 9 08:30:52 1999 From: o.kirstein at fz-juelich.de (Dr. Oliver Kirstein) Date: Wed, 09 Jun 1999 08:30:52 +0200 Subject: Questions Message-ID: <375E0A1C.FC6F3781@fz-juelich.de> Hi Kristian! I'm still working on the PST component but I have now another question concerning the handling of the neutron coordinates in McStas: I try to build a large focusing monochromator (horizontal & vertical) for a backscattering device (this is simply a three axis spectrometer where the third axis points in the opposite direction of the first defelector crystal. I attached two examples focus69.instr (somthing like a 3-axis spec) and focus90.instr (backscattering type) and you will se the difference with mcdisplay. What make me wonder is the fact that the number of neutrons which reach the final energy monitor reduces drastically in the case of focus90 although the reflectivity due to the 90-deg-geometry should be better. Is it possible that the neutron passes the first deflector crystal (called mono1) twice in contrast to the focus69-type and that the neutrons detected are those which have been transmitted? (I started both programs with focusXX --ncount 1.e5 --dir=simXX --seed=1) Thank you for your help. Regards, Oliver -- =============================================================================== _/_/_/ _/_/_/_/_/ _/_/_/_/_/ Oliver Kirstein _/ _/ _/ IFF, Neutronenstreuung _/ _/_/_/ _/_/_/ Forschungszentrum Juelich _/ _/ _/ D-52425 Juelich _/ _/ _/ Tel.: +49 (0)2461 614532 _/_/_/ _/ _/ Fax.: +49 (0)2461 612610 http://iff034.iff.kfa-juelich.de/Oliver/ =============================================================================== -------------- next part -------------- DEFINE INSTRUMENT focus() DECLARE %{ double TTM, TTM2; double EN, DEN, LMIN, LMAX; double ZMINM, ZMAXM, YMINM, YMAXM, MOSAICHM, MOSAICVM, R0M, QMM; double OMM, MONO_Y, MONO_Z; double OMM2, Q2; // #include "mono_var" double XPOSM_0 , YPOSM_0 , ZPOSM_0 ; double XROTM_0 , YROTM_0 , ZROTM_0 ; double XPOSM_1 , YPOSM_1 , ZPOSM_1 ; double XROTM_1 , YROTM_1 , ZROTM_1 ; double XPOSM_2 , YPOSM_2 , ZPOSM_2 ; double XROTM_2 , YROTM_2 , ZROTM_2 ; double XPOSM_3 , YPOSM_3 , ZPOSM_3 ; double XROTM_3 , YROTM_3 , ZROTM_3 ; double XPOSM_4 , YPOSM_4 , ZPOSM_4 ; double XROTM_4 , YROTM_4 , ZROTM_4 ; double XPOSM_5 , YPOSM_5 , ZPOSM_5 ; double XROTM_5 , YROTM_5 , ZROTM_5 ; double XPOSM_6 , YPOSM_6 , ZPOSM_6 ; double XROTM_6 , YROTM_6 , ZROTM_6 ; double XPOSM_7 , YPOSM_7 , ZPOSM_7 ; double XROTM_7 , YROTM_7 , ZROTM_7 ; double XPOSM_8 , YPOSM_8 , ZPOSM_8 ; double XROTM_8 , YROTM_8 , ZROTM_8 ; double XPOSM_9 , YPOSM_9 , ZPOSM_9 ; double XROTM_9 , YROTM_9 , ZROTM_9 ; double XPOSM_10 , YPOSM_10 , ZPOSM_10 ; double XROTM_10 , YROTM_10 , ZROTM_10 ; double XPOSM_11 , YPOSM_11 , ZPOSM_11 ; double XROTM_11 , YROTM_11 , ZROTM_11 ; double XPOSM_12 , YPOSM_12 , ZPOSM_12 ; double XROTM_12 , YROTM_12 , ZROTM_12 ; double XPOSM_13 , YPOSM_13 , ZPOSM_13 ; double XROTM_13 , YROTM_13 , ZROTM_13 ; double XPOSM_14 , YPOSM_14 , ZPOSM_14 ; double XROTM_14 , YROTM_14 , ZROTM_14 ; double XPOSM_15 , YPOSM_15 , ZPOSM_15 ; double XROTM_15 , YROTM_15 , ZROTM_15 ; double XPOSM_16 , YPOSM_16 , ZPOSM_16 ; double XROTM_16 , YROTM_16 , ZROTM_16 ; double XPOSM_17 , YPOSM_17 , ZPOSM_17 ; double XROTM_17 , YROTM_17 , ZROTM_17 ; double XPOSM_18 , YPOSM_18 , ZPOSM_18 ; double XROTM_18 , YROTM_18 , ZROTM_18 ; double XPOSM_19 , YPOSM_19 , ZPOSM_19 ; double XROTM_19 , YROTM_19 , ZROTM_19 ; double XPOSM_20 , YPOSM_20 , ZPOSM_20 ; double XROTM_20 , YROTM_20 , ZROTM_20 ; %} INITIALIZE %{ //#include "mono_pos" XPOSM_0 = -0.24913305, XROTM_0 = -2.38732409; YPOSM_0 = -0.0833092257, YROTM_0 = 60.0380249; ZPOSM_0 = 2.18267322, ZROTM_0 = 0.; XPOSM_1 = -0.24934946, XROTM_1 = 0.; YPOSM_1 = 0., YROTM_1 = 60.0380249; ZPOSM_1 = 2.18439531, ZROTM_1 = 0.; XPOSM_2 = -0.24913305, XROTM_2 = 2.38732409; YPOSM_2 = 0.0833092257, YROTM_2 = 60.0380249; ZPOSM_2 = 2.18267322, ZROTM_2 = 0.; XPOSM_3 = -0.166329354, XROTM_3 = -2.38732409; YPOSM_3 = -0.0833092257, YROTM_3 = 62.4253502; ZPOSM_3 = 2.19132972, ZROTM_3 = 0.; XPOSM_4 = -0.166473836, XROTM_4 = 0.; YPOSM_4 = 0., YROTM_4 = 62.4253502; ZPOSM_4 = 2.19305944, ZROTM_4 = 0.; XPOSM_5 = -0.166329354, XROTM_5 = 2.38732409; YPOSM_5 = 0.0833092257, YROTM_5 = 62.4253502; ZPOSM_5 = 2.19132972, ZROTM_5 = 0.; XPOSM_6 = -0.0832369179, XROTM_6 = -2.38732409; YPOSM_6 = -0.0833092257, YROTM_6 = 64.8126755; ZPOSM_6 = 2.19652987, ZROTM_6 = 0.; XPOSM_7 = -0.0833092257, XROTM_7 = 0.; YPOSM_7 = 0., YROTM_7 = 64.8126755; ZPOSM_7 = 2.19826412, ZROTM_7 = 0.; XPOSM_8 = -0.0832369179, XROTM_8 = 2.38732409; YPOSM_8 = 0.0833092257, YROTM_8 = 64.8126755; ZPOSM_8 = 2.19652987, ZROTM_8 = 0.; XPOSM_9 = 0., XROTM_9 = -2.38732409; YPOSM_9 = -0.0833092257, YROTM_9 = 67.1999969; ZPOSM_9 = 2.19826412, ZROTM_9 = 0.; XPOSM_10 = 0., XROTM_10 = 0.; YPOSM_10 = 0., YROTM_10 = 67.1999969; ZPOSM_10 = 2.20000005, ZROTM_10 = 0.; XPOSM_11 = 0., XROTM_11 = 2.38732409; YPOSM_11 = 0.0833092257, YROTM_11 = 67.1999969; ZPOSM_11 = 2.19826412, ZROTM_11 = 0.; XPOSM_12 = 0.0832369179, XROTM_12 = -2.38732409; YPOSM_12 = -0.0833092257, YROTM_12 = 69.5873184; ZPOSM_12 = 2.19652987, ZROTM_12 = 0.; XPOSM_13 = 0.0833092257, XROTM_13 = 0.; YPOSM_13 = 0., YROTM_13 = 69.5873184; ZPOSM_13 = 2.19826412, ZROTM_13 = 0.; XPOSM_14 = 0.0832369179, XROTM_14 = 2.38732409; YPOSM_14 = 0.0833092257, YROTM_14 = 69.5873184; ZPOSM_14 = 2.19652987, ZROTM_14 = 0.; XPOSM_15 = 0.166329354, XROTM_15 = -2.38732409; YPOSM_15 = -0.0833092257, YROTM_15 = 71.9746475; ZPOSM_15 = 2.19132972, ZROTM_15 = 0.; XPOSM_16 = 0.166473836, XROTM_16 = 0.; YPOSM_16 = 0., YROTM_16 = 71.9746475; ZPOSM_16 = 2.19305944, ZROTM_16 = 0.; XPOSM_17 = 0.166329354, XROTM_17 = 2.38732409; YPOSM_17 = 0.0833092257, YROTM_17 = 71.9746475; ZPOSM_17 = 2.19132972, ZROTM_17 = 0.; XPOSM_18 = 0.24913305, XROTM_18 = -2.38732409; YPOSM_18 = -0.0833092257, YROTM_18 = 74.361969; ZPOSM_18 = 2.18267322, ZROTM_18 = 0.; XPOSM_19 = 0.24934946, XROTM_19 = 0.; YPOSM_19 = 0., YROTM_19 = 74.361969; ZPOSM_19 = 2.18439531, ZROTM_19 = 0.; XPOSM_20 = 0.24913305, XROTM_20 = 2.38732409; YPOSM_20 = 0.0833092257, YROTM_20 = 74.361969; ZPOSM_20 = 2.18267322, ZROTM_20 = 0.; OMM = 69.2; OMM2 = 69.2; Q2 = 1.8733688; MONO_Y = 0.0417; MONO_Z = 0.0417; ZMINM = -MONO_Z; ZMAXM = MONO_Z; YMINM = -MONO_Y; YMAXM = MONO_Y; MOSAICHM = 60.0; MOSAICVM = 60.0; R0M = 1.0; QMM = Q2; TTM = 2.0 * OMM; TTM2 = 2.0 * OMM2; EN = 2.08; DEN = 0.24; LMIN = sqrt(81.81 / (EN + DEN) ); LMAX = sqrt(81.81 / (EN - DEN) ); %} TRACE COMPONENT a1 = Arm() AT (0,0,0) ABSOLUTE COMPONENT source = Source_flux(radius = 0.01, dist = 2.0, xw=0.05, yh=0.05, E0=EN, dE=DEN, flux=8.0e+14) AT (0, 0, 0) RELATIVE a1 ROTATED (0,0,0) RELATIVE a1 COMPONENT Slit1 = Slit (xmin=-0.03, xmax= 0.03, ymin=-0.03, ymax=0.03) AT (0, 0, 0.5) RELATIVE a1 ROTATED (0,0,0) RELATIVE a1 COMPONENT E0 = E_monitor (xmin=-0.03, xmax=0.03, ymin=-0.03, ymax=0.03, Emin=1.8, Emax=2.3, nchan=100, filename="e0") AT (0,0,2.0) RELATIVE a1 ROTATED (0,0,0) RELATIVE a1 COMPONENT mono1 = Monochromator (zmin=-0.03, zmax=0.03, ymin=-0.015, ymax=0.015 , mosaich=240.0, mosaicv=60.0, r0=1.0, Q=1.8733688) AT (0,0,2.0) RELATIVE a1 ROTATED (0,OMM,0) RELATIVE a1 COMPONENT a2 = Arm() AT (0,0,0) RELATIVE mono1 ROTATED (0,TTM,0) RELATIVE a1 COMPONENT E1 = E_monitor (xmin=-0.03, xmax=0.03, ymin=-0.03, ymax=0.03, Emin=1.8, Emax=2.3, nchan=100, filename="e1") AT (0,0,0.1) RELATIVE a2 ROTATED (0,0,0) RELATIVE a2 //%include "mono_def" /***************************************************/ /*** begin of monochromator definition ***/ /*** ***/ component Doppler_0 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_0 ,YPOSM_0 ,ZPOSM_0 ) relative a2 rotated (XROTM_0 ,YROTM_0 ,ZROTM_0 ) relative a2 component Doppler_1 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_1 ,YPOSM_1 ,ZPOSM_1 ) relative a2 rotated (XROTM_1 ,YROTM_1 ,ZROTM_1 ) relative a2 component Doppler_2 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_2 ,YPOSM_2 ,ZPOSM_2 ) relative a2 rotated (XROTM_2 ,YROTM_2 ,ZROTM_2 ) relative a2 component Doppler_3 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_3 ,YPOSM_3 ,ZPOSM_3 ) relative a2 rotated (XROTM_3 ,YROTM_3 ,ZROTM_3 ) relative a2 component Doppler_4 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_4 ,YPOSM_4 ,ZPOSM_4 ) relative a2 rotated (XROTM_4 ,YROTM_4 ,ZROTM_4 ) relative a2 component Doppler_5 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_5 ,YPOSM_5 ,ZPOSM_5 ) relative a2 rotated (XROTM_5 ,YROTM_5 ,ZROTM_5 ) relative a2 component Doppler_6 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_6 ,YPOSM_6 ,ZPOSM_6 ) relative a2 rotated (XROTM_6 ,YROTM_6 ,ZROTM_6 ) relative a2 component Doppler_7 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_7 ,YPOSM_7 ,ZPOSM_7 ) relative a2 rotated (XROTM_7 ,YROTM_7 ,ZROTM_7 ) relative a2 component Doppler_8 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_8 ,YPOSM_8 ,ZPOSM_8 ) relative a2 rotated (XROTM_8 ,YROTM_8 ,ZROTM_8 ) relative a2 component Doppler_9 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_9 ,YPOSM_9 ,ZPOSM_9 ) relative a2 rotated (XROTM_9 ,YROTM_9 ,ZROTM_9 ) relative a2 component Doppler_10 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_10 ,YPOSM_10 ,ZPOSM_10 ) relative a2 rotated (XROTM_10 ,YROTM_10 ,ZROTM_10 ) relative a2 component Doppler_11 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_11 ,YPOSM_11 ,ZPOSM_11 ) relative a2 rotated (XROTM_11 ,YROTM_11 ,ZROTM_11 ) relative a2 component Doppler_12 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_12 ,YPOSM_12 ,ZPOSM_12 ) relative a2 rotated (XROTM_12 ,YROTM_12 ,ZROTM_12 ) relative a2 component Doppler_13 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_13 ,YPOSM_13 ,ZPOSM_13 ) relative a2 rotated (XROTM_13 ,YROTM_13 ,ZROTM_13 ) relative a2 component Doppler_14 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_14 ,YPOSM_14 ,ZPOSM_14 ) relative a2 rotated (XROTM_14 ,YROTM_14 ,ZROTM_14 ) relative a2 component Doppler_15 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_15 ,YPOSM_15 ,ZPOSM_15 ) relative a2 rotated (XROTM_15 ,YROTM_15 ,ZROTM_15 ) relative a2 component Doppler_16 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_16 ,YPOSM_16 ,ZPOSM_16 ) relative a2 rotated (XROTM_16 ,YROTM_16 ,ZROTM_16 ) relative a2 component Doppler_17 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_17 ,YPOSM_17 ,ZPOSM_17 ) relative a2 rotated (XROTM_17 ,YROTM_17 ,ZROTM_17 ) relative a2 component Doppler_18 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_18 ,YPOSM_18 ,ZPOSM_18 ) relative a2 rotated (XROTM_18 ,YROTM_18 ,ZROTM_18 ) relative a2 component Doppler_19 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_19 ,YPOSM_19 ,ZPOSM_19 ) relative a2 rotated (XROTM_19 ,YROTM_19 ,ZROTM_19 ) relative a2 component Doppler_20 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_20 ,YPOSM_20 ,ZPOSM_20 ) relative a2 rotated (XROTM_20 ,YROTM_20 ,ZROTM_20 ) relative a2 component a3 = Arm() at (0,0,0) relative Doppler_10 rotated (0,TTM2,0) relative a2 /*** end of monochromator definition ***/ /*** ***/ /***************************************************/ COMPONENT E2 = E_monitor (xmin=-0.015, xmax=0.015, ymin=-0.015, ymax=0.015, Emin=1.8, Emax=2.3, nchan=100, filename="e2") AT (0,0,2.4) RELATIVE a3 ROTATED (0,0,0) RELATIVE a3 END -------------- next part -------------- DEFINE INSTRUMENT focus90() DECLARE %{ double TTM, TTM2; double EN, DEN, LMIN, LMAX; double ZMINM, ZMAXM, YMINM, YMAXM, MOSAICHM, MOSAICVM, R0M, QMM; double OMM, MONO_Y, MONO_Z; double OMM2, Q2; //#include "mono_var" double XPOSM_0 , YPOSM_0 , ZPOSM_0 ; double XROTM_0 , YROTM_0 , ZROTM_0 ; double XPOSM_1 , YPOSM_1 , ZPOSM_1 ; double XROTM_1 , YROTM_1 , ZROTM_1 ; double XPOSM_2 , YPOSM_2 , ZPOSM_2 ; double XROTM_2 , YROTM_2 , ZROTM_2 ; double XPOSM_3 , YPOSM_3 , ZPOSM_3 ; double XROTM_3 , YROTM_3 , ZROTM_3 ; double XPOSM_4 , YPOSM_4 , ZPOSM_4 ; double XROTM_4 , YROTM_4 , ZROTM_4 ; double XPOSM_5 , YPOSM_5 , ZPOSM_5 ; double XROTM_5 , YROTM_5 , ZROTM_5 ; double XPOSM_6 , YPOSM_6 , ZPOSM_6 ; double XROTM_6 , YROTM_6 , ZROTM_6 ; double XPOSM_7 , YPOSM_7 , ZPOSM_7 ; double XROTM_7 , YROTM_7 , ZROTM_7 ; double XPOSM_8 , YPOSM_8 , ZPOSM_8 ; double XROTM_8 , YROTM_8 , ZROTM_8 ; double XPOSM_9 , YPOSM_9 , ZPOSM_9 ; double XROTM_9 , YROTM_9 , ZROTM_9 ; double XPOSM_10 , YPOSM_10 , ZPOSM_10 ; double XROTM_10 , YROTM_10 , ZROTM_10 ; double XPOSM_11 , YPOSM_11 , ZPOSM_11 ; double XROTM_11 , YROTM_11 , ZROTM_11 ; double XPOSM_12 , YPOSM_12 , ZPOSM_12 ; double XROTM_12 , YROTM_12 , ZROTM_12 ; double XPOSM_13 , YPOSM_13 , ZPOSM_13 ; double XROTM_13 , YROTM_13 , ZROTM_13 ; double XPOSM_14 , YPOSM_14 , ZPOSM_14 ; double XROTM_14 , YROTM_14 , ZROTM_14 ; double XPOSM_15 , YPOSM_15 , ZPOSM_15 ; double XROTM_15 , YROTM_15 , ZROTM_15 ; double XPOSM_16 , YPOSM_16 , ZPOSM_16 ; double XROTM_16 , YROTM_16 , ZROTM_16 ; double XPOSM_17 , YPOSM_17 , ZPOSM_17 ; double XROTM_17 , YROTM_17 , ZROTM_17 ; double XPOSM_18 , YPOSM_18 , ZPOSM_18 ; double XROTM_18 , YROTM_18 , ZROTM_18 ; double XPOSM_19 , YPOSM_19 , ZPOSM_19 ; double XROTM_19 , YROTM_19 , ZROTM_19 ; double XPOSM_20 , YPOSM_20 , ZPOSM_20 ; double XROTM_20 , YROTM_20 , ZROTM_20 ; %} INITIALIZE %{ //#include "mono_pos" XPOSM_0 = -0.24913305, XROTM_0 = -2.38732409; YPOSM_0 = -0.0833092257, YROTM_0 = 82.838028; ZPOSM_0 = 2.18267322, ZROTM_0 = 0.; XPOSM_1 = -0.24934946, XROTM_1 = 0.; YPOSM_1 = 0., YROTM_1 = 82.838028; ZPOSM_1 = 2.18439531, ZROTM_1 = 0.; XPOSM_2 = -0.24913305, XROTM_2 = 2.38732409; YPOSM_2 = 0.0833092257, YROTM_2 = 82.838028; ZPOSM_2 = 2.18267322, ZROTM_2 = 0.; XPOSM_3 = -0.166329354, XROTM_3 = -2.38732409; YPOSM_3 = -0.0833092257, YROTM_3 = 85.2253494; ZPOSM_3 = 2.19132972, ZROTM_3 = 0.; XPOSM_4 = -0.166473836, XROTM_4 = 0.; YPOSM_4 = 0., YROTM_4 = 85.2253494; ZPOSM_4 = 2.19305944, ZROTM_4 = 0.; XPOSM_5 = -0.166329354, XROTM_5 = 2.38732409; YPOSM_5 = 0.0833092257, YROTM_5 = 85.2253494; ZPOSM_5 = 2.19132972, ZROTM_5 = 0.; XPOSM_6 = -0.0832369179, XROTM_6 = -2.38732409; YPOSM_6 = -0.0833092257, YROTM_6 = 87.6126785; ZPOSM_6 = 2.19652987, ZROTM_6 = 0.; XPOSM_7 = -0.0833092257, XROTM_7 = 0.; YPOSM_7 = 0., YROTM_7 = 87.6126785; ZPOSM_7 = 2.19826412, ZROTM_7 = 0.; XPOSM_8 = -0.0832369179, XROTM_8 = 2.38732409; YPOSM_8 = 0.0833092257, YROTM_8 = 87.6126785; ZPOSM_8 = 2.19652987, ZROTM_8 = 0.; XPOSM_9 = 0., XROTM_9 = -2.38732409; YPOSM_9 = -0.0833092257, YROTM_9 = 90.; ZPOSM_9 = 2.19826412, ZROTM_9 = 0.; XPOSM_10 = 0., XROTM_10 = 0.; YPOSM_10 = 0., YROTM_10 = 90.; ZPOSM_10 = 2.20000005, ZROTM_10 = 0.; XPOSM_11 = 0., XROTM_11 = 2.38732409; YPOSM_11 = 0.0833092257, YROTM_11 = 90.; ZPOSM_11 = 2.19826412, ZROTM_11 = 0.; XPOSM_12 = 0.0832369179, XROTM_12 = -2.38732409; YPOSM_12 = -0.0833092257, YROTM_12 = 92.3873215; ZPOSM_12 = 2.19652987, ZROTM_12 = 0.; XPOSM_13 = 0.0833092257, XROTM_13 = 0.; YPOSM_13 = 0., YROTM_13 = 92.3873215; ZPOSM_13 = 2.19826412, ZROTM_13 = 0.; XPOSM_14 = 0.0832369179, XROTM_14 = 2.38732409; YPOSM_14 = 0.0833092257, YROTM_14 = 92.3873215; ZPOSM_14 = 2.19652987, ZROTM_14 = 0.; XPOSM_15 = 0.166329354, XROTM_15 = -2.38732409; YPOSM_15 = -0.0833092257, YROTM_15 = 94.7746506; ZPOSM_15 = 2.19132972, ZROTM_15 = 0.; XPOSM_16 = 0.166473836, XROTM_16 = 0.; YPOSM_16 = 0., YROTM_16 = 94.7746506; ZPOSM_16 = 2.19305944, ZROTM_16 = 0.; XPOSM_17 = 0.166329354, XROTM_17 = 2.38732409; YPOSM_17 = 0.0833092257, YROTM_17 = 94.7746506; ZPOSM_17 = 2.19132972, ZROTM_17 = 0.; XPOSM_18 = 0.24913305, XROTM_18 = -2.38732409; YPOSM_18 = -0.0833092257, YROTM_18 = 97.161972; ZPOSM_18 = 2.18267322, ZROTM_18 = 0.; XPOSM_19 = 0.24934946, XROTM_19 = 0.; YPOSM_19 = 0., YROTM_19 = 97.161972; ZPOSM_19 = 2.18439531, ZROTM_19 = 0.; XPOSM_20 = 0.24913305, XROTM_20 = 2.38732409; YPOSM_20 = 0.0833092257, YROTM_20 = 97.161972; ZPOSM_20 = 2.18267322, ZROTM_20 = 0.; OMM = 69.2; OMM2 = 90.0; Q2 = 2.004; MONO_Y = 0.0417; MONO_Z = 0.0417; ZMINM = -MONO_Z; ZMAXM = MONO_Z; YMINM = -MONO_Y; YMAXM = MONO_Y; MOSAICHM = 60.0; MOSAICVM = 60.0; R0M = 1.0; QMM = Q2; TTM = 2.0 * OMM; TTM2 = 2.0 * OMM2; EN = 2.08; DEN = 0.24; LMIN = sqrt(81.81 / (EN + DEN) ); LMAX = sqrt(81.81 / (EN - DEN) ); %} TRACE COMPONENT a1 = Arm() AT (0,0,0) ABSOLUTE COMPONENT source = Source_flux(radius = 0.01, dist = 2.0, xw=0.05, yh=0.05, E0=EN, dE=DEN, flux=8.0e+14) AT (0, 0, 0) RELATIVE a1 ROTATED (0,0,0) RELATIVE a1 COMPONENT Slit1 = Slit (xmin=-0.03, xmax= 0.03, ymin=-0.03, ymax=0.03) AT (0, 0, 0.5) RELATIVE a1 ROTATED (0,0,0) RELATIVE a1 COMPONENT E0 = E_monitor (xmin=-0.03, xmax=0.03, ymin=-0.03, ymax=0.03, Emin=1.8, Emax=2.3, nchan=100, filename="e0") AT (0,0,2.0) RELATIVE a1 ROTATED (0,0,0) RELATIVE a1 COMPONENT mono1 = Monochromator (zmin=-0.03, zmax=0.03, ymin=-0.015, ymax=0.015 , mosaich=240.0, mosaicv=60.0, r0=1.0, Q=1.8733688) AT (0,0,2.0) RELATIVE a1 ROTATED (0,OMM,0) RELATIVE a1 COMPONENT a2 = Arm() AT (0,0,0) RELATIVE mono1 ROTATED (0,TTM,0) RELATIVE a1 COMPONENT E1 = E_monitor (xmin=-0.03, xmax=0.03, ymin=-0.03, ymax=0.03, Emin=1.8, Emax=2.3, nchan=100, filename="e1") AT (0,0,0.1) RELATIVE a2 ROTATED (0,0,0) RELATIVE a2 //%include "mono_def" /***************************************************/ /*** begin of monochromator definition ***/ /*** ***/ component Doppler_0 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_0 ,YPOSM_0 ,ZPOSM_0 ) relative a2 rotated (XROTM_0 ,YROTM_0 ,ZROTM_0 ) relative a2 component Doppler_1 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_1 ,YPOSM_1 ,ZPOSM_1 ) relative a2 rotated (XROTM_1 ,YROTM_1 ,ZROTM_1 ) relative a2 component Doppler_2 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_2 ,YPOSM_2 ,ZPOSM_2 ) relative a2 rotated (XROTM_2 ,YROTM_2 ,ZROTM_2 ) relative a2 component Doppler_3 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_3 ,YPOSM_3 ,ZPOSM_3 ) relative a2 rotated (XROTM_3 ,YROTM_3 ,ZROTM_3 ) relative a2 component Doppler_4 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_4 ,YPOSM_4 ,ZPOSM_4 ) relative a2 rotated (XROTM_4 ,YROTM_4 ,ZROTM_4 ) relative a2 component Doppler_5 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_5 ,YPOSM_5 ,ZPOSM_5 ) relative a2 rotated (XROTM_5 ,YROTM_5 ,ZROTM_5 ) relative a2 component Doppler_6 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_6 ,YPOSM_6 ,ZPOSM_6 ) relative a2 rotated (XROTM_6 ,YROTM_6 ,ZROTM_6 ) relative a2 component Doppler_7 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_7 ,YPOSM_7 ,ZPOSM_7 ) relative a2 rotated (XROTM_7 ,YROTM_7 ,ZROTM_7 ) relative a2 component Doppler_8 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_8 ,YPOSM_8 ,ZPOSM_8 ) relative a2 rotated (XROTM_8 ,YROTM_8 ,ZROTM_8 ) relative a2 component Doppler_9 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_9 ,YPOSM_9 ,ZPOSM_9 ) relative a2 rotated (XROTM_9 ,YROTM_9 ,ZROTM_9 ) relative a2 component Doppler_10 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_10 ,YPOSM_10 ,ZPOSM_10 ) relative a2 rotated (XROTM_10 ,YROTM_10 ,ZROTM_10 ) relative a2 component Doppler_11 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_11 ,YPOSM_11 ,ZPOSM_11 ) relative a2 rotated (XROTM_11 ,YROTM_11 ,ZROTM_11 ) relative a2 component Doppler_12 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_12 ,YPOSM_12 ,ZPOSM_12 ) relative a2 rotated (XROTM_12 ,YROTM_12 ,ZROTM_12 ) relative a2 component Doppler_13 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_13 ,YPOSM_13 ,ZPOSM_13 ) relative a2 rotated (XROTM_13 ,YROTM_13 ,ZROTM_13 ) relative a2 component Doppler_14 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_14 ,YPOSM_14 ,ZPOSM_14 ) relative a2 rotated (XROTM_14 ,YROTM_14 ,ZROTM_14 ) relative a2 component Doppler_15 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_15 ,YPOSM_15 ,ZPOSM_15 ) relative a2 rotated (XROTM_15 ,YROTM_15 ,ZROTM_15 ) relative a2 component Doppler_16 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_16 ,YPOSM_16 ,ZPOSM_16 ) relative a2 rotated (XROTM_16 ,YROTM_16 ,ZROTM_16 ) relative a2 component Doppler_17 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_17 ,YPOSM_17 ,ZPOSM_17 ) relative a2 rotated (XROTM_17 ,YROTM_17 ,ZROTM_17 ) relative a2 component Doppler_18 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_18 ,YPOSM_18 ,ZPOSM_18 ) relative a2 rotated (XROTM_18 ,YROTM_18 ,ZROTM_18 ) relative a2 component Doppler_19 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_19 ,YPOSM_19 ,ZPOSM_19 ) relative a2 rotated (XROTM_19 ,YROTM_19 ,ZROTM_19 ) relative a2 component Doppler_20 = Monochromator(zmin=ZMINM, zmax=ZMAXM, ymin=YMINM, ymax=YMAXM, mosaich=MOSAICHM, mosaicv=MOSAICVM, r0=R0M, Q=QMM) at (XPOSM_20 ,YPOSM_20 ,ZPOSM_20 ) relative a2 rotated (XROTM_20 ,YROTM_20 ,ZROTM_20 ) relative a2 component a3 = Arm() at (0,0,0) relative Doppler_10 rotated (0,TTM2,0) relative a2 /*** end of monochromator definition ***/ /*** ***/ /***************************************************/ COMPONENT E2 = E_monitor (xmin=-0.015, xmax=0.015, ymin=-0.015, ymax=0.015, Emin=1.8, Emax=2.3, nchan=100, filename="e2") AT (0,0,2.4) RELATIVE a3 ROTATED (0,0,0) RELATIVE a3 END -------------- next part -------------- A non-text attachment was scrubbed... Name: vcard.vcf Type: text/x-vcard Size: 480 bytes Desc: Card for Dr. Oliver Kirstein URL: From kristian.nielsen at risoe.dk Wed Jun 9 11:22:46 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: 9 Jun 1999 11:22:46 +0200 Subject: Questions In-Reply-To: <375E0A1C.FC6F3781@fz-juelich.de> (o.kirstein@fz-juelich.de) Message-ID: > I'm still working on the PST component but I have now another question > concerning the handling of the neutron coordinates in McStas: I try to > build a large focusing monochromator (horizontal & vertical) for a > backscattering device (this is simply a three axis spectrometer where > the third axis points in the opposite direction of the first defelector > crystal. I attached two examples focus69.instr (somthing like a 3-axis > spec) and focus90.instr (backscattering type) and you will se the > difference with mcdisplay. Ok, I tried your examples and did a little experimentation. First of all, there is a bug in the Monochromator component. When the neutron wavelength is longer than the maximum for which bragg scattering is possible, an invalid operation (asin(x) with x > 1) is performed. This seems to happen in your focus90 instrument. I have a corrected component, and I just made it available from the "known bugs" page: http://neutron.risoe.dk/mcstas/known-bugs.html But I am not sure if this problem affects your results. > What make me wonder is the fact that the number of neutrons which reach > the final energy monitor reduces drastically in the case of focus90 > although the reflectivity due to the 90-deg-geometry should be better. The number of neutrons in the detector generally carries no physical significance in McStas. If you look instead at the intensities, they are 1.8e6 for focus69 and 1.3e7 (8 times bigger) for focus90, so this is as you expected, I think? The decrease in the number of neutrons is of course a problem for the efficiency of the simulation. The reason is that you have a very broad energy band emission from the source compared to the very small energy band accepted by the large monochromator; most of the neutrons have an energy that makes it impossible to ever reach the final detector. (I actually have a new source component that automatically adapts itself to situations like this and picks more neutrons with useful energies, but this is still a bit experimental. Meanwhile, if you match the energy range in the source to the range accepted by the large monochromator, you should get better statistics.) Hm, actually the accepted energy range in focus69 is much bigger than in the backscattering focus90 mode, so I would actually have expected the intensity in backscattering to be lower than at 69 degrees. But I guess that there is a simple explanation for this? > Is it possible that the neutron passes the first deflector crystal > (called mono1) twice in contrast to the focus69-type and that the > neutrons detected are those which have been transmitted? No. A component can interact at most once with each neutron, regardless of its physical position. I hope this answers your questions; otherwise please ask again. Things are currently a bit quiet with regards to McStas, so I am happy to know that it is still being used. - Kristian. From o.kirstein at fz-juelich.de Wed Jun 9 13:35:05 1999 From: o.kirstein at fz-juelich.de (Dr. Oliver Kirstein) Date: Wed, 09 Jun 1999 13:35:05 +0200 Subject: Questions References: <01JC70Q3ASRK94KRWH@risoe.dk> Message-ID: <375E5169.21D45C73@fz-juelich.de> First off all: Thanks for the reply (you're right with the intensities, I was a little bit too fast). Another point concerning the monochromator component: The analyser component of focus90/69 consists of perfect crystals which means that they normally have no mosaic spread (mosaich=mosaicv=0). While the calculation of the mc-probability the p-value is calculated via tmp3 (prop. to 1/mosaich) so p should be set to 1 in this case (I hope that I'm right). Bye, Oliver -- =============================================================================== _/_/_/ _/_/_/_/_/ _/_/_/_/_/ Dr. Oliver Kirstein _/ _/ _/ IFF, Neutronenstreuung _/ _/_/_/ _/_/_/ Forschungszentrum Juelich _/ _/ _/ D-52425 Juelich _/ _/ _/ Tel.: +49 (0)2461 614532 _/_/_/ _/ _/ Fax.: +49 (0)2461 612610 http://iff034.iff.kfa-juelich.de/Oliver/ =============================================================================== -------------- next part -------------- A non-text attachment was scrubbed... Name: vcard.vcf Type: text/x-vcard Size: 480 bytes Desc: Card for Dr. Oliver Kirstein URL: From Bertrand.Roessli at psi.ch Wed Jun 9 15:26:54 1999 From: Bertrand.Roessli at psi.ch (B.Roessli) Date: Wed, 09 Jun 1999 14:26:54 +0100 Subject: Mcstas 1.1 Message-ID: <375E6B98.6AAEB7AC@psi.ch> Dear Dr. Nielsen, We installed the new Mcstas 1.1 version at PSI. We have the following problems 1. After compiling the program for the linup-1.instr example, it turns out that the globale double variables "mns1" and "mns2' are not declared. 2. using cc with DU 4.OD gives a core dump. The prgram works with gcc. Thank You for helping us, regards, B. Roessli From kristian.nielsen at risoe.dk Wed Jun 9 15:22:18 1999 From: kristian.nielsen at risoe.dk (Kristian Nielsen) Date: 9 Jun 1999 15:22:18 +0200 Subject: Mcstas 1.1 In-Reply-To: <375E6B98.6AAEB7AC@psi.ch> (Bertrand.Roessli@psi.ch) Message-ID: > We installed the new Mcstas 1.1 version at PSI. > We have the following problems > 1. After compiling the program for the linup-1.instr example, it turns > out that the globale double variables > "mns1" and "mns2' are not declared. Hm, I do not understand this. Maybe you could give more details? Could you send me a copy of the "linup-1.c" file produced by McStas with this problem? I tried to log in on my account on your host "lnsa15.psi.ch" and run the linup-1 example. The version of McStas installed on this machine appears to be v1.0, and I could compile and run the linup-1 example without problems. I also had no problems with my private installation of McStas v1.1 on this machine. > 2. using cc with DU 4.OD gives a core dump. The prgram works with gcc. Yes, this is a known problem/feature of the Digital Unix C compiler. I just added an explanation of this problem to a new "FAQ" with this address: http://neutron.risoe.dk/mcstas/FAQ.html Basically, it is necessary to pass the "-std1" option to the Digital C compiler to get ANSI C semantics, like this: cc -std1 -o linup-1 linup-1.c -lm > Thank You for helping us, You are most welcome! If you have any further problems, feel free to ask again. - Kristian. From Bertrand.Roessli at psi.ch Thu Jun 10 13:53:56 1999 From: Bertrand.Roessli at psi.ch (Bertrand.Roessli at psi.ch) Date: Thu, 10 Jun 1999 13:53:56 +0200 Subject: Mcstas 1.1 Message-ID: <99061013535629@psicl0.psi.ch> Dear Kristian, please find below the linup-1.c file which has got undefined variables. regards, B. Roessli /* Automatically generated file. Do not edit. */ #define MC_USE_DEFAULT_MAIN #define MC_EMBEDDED_RUNTIME #line 1 "mcstas-r.h" /******************************************************************************* * Runtime system for McStas. * * Project: Monte Carlo Simulation of Triple Axis Spectrometers * File name: mcstas-r.h * * Author: K.N. Aug 29, 1997 * * $Id: mcstas-r.h,v 1.20 1998/10/09 07:53:48 kn Exp kn $ * * $Log: mcstas-r.h,v $ * Revision 1.20 1998/10/09 07:53:48 kn * Added some unit conversion constants. * * Revision 1.19 1998/10/02 08:38:36 kn * Added DETECTOR_OUT support. * Fixed header comment. * * Revision 1.18 1998/10/01 08:12:42 kn * Support for embedding the file in the output from McStas. * Added mcstas_main() function. * Added support for command line arguments. * * Revision 1.17 1998/09/24 13:01:39 kn * Minor conversion factor additions. * * Revision 1.16 1998/09/23 13:52:08 kn * Added conversion factors. * McStas now uses its own random() implementation (unless * USE_SYSTEM_RANDOM is defined). * * Revision 1.15 1998/05/19 07:59:45 kn * Hack to make random number generation work with HP's CC C compiler. * * Revision 1.14 1998/04/17 11:50:31 kn * Added sphere_intersect. * * Revision 1.13 1998/04/17 10:53:08 kn * Added randvec_target_sphere. * * Revision 1.12 1998/03/25 07:23:24 kn * Fixed RAND_MAX #define for HPUX. * * Revision 1.11 1998/03/24 13:59:26 lefmann * Added #define for RAND_MAX, needed on HPUX. * * Revision 1.10 1998/03/24 13:24:40 lefmann * Added HBAR, MNEUTRON. * * Revision 1.9 1998/03/24 07:42:35 kn * Added definition for PI. * * Revision 1.8 1998/03/24 07:36:20 kn * Make ABSORB macro work better in control structures. * Add test_printf(). * Add rand01(), randpm1(), rand0max(), and randminmax(). * Add PROP_X0, PROP_Y0, PROP_DT, vec_prod(), scalar_prod(), NORM(), and * rotate(). * Fix typos. * * Revision 1.7 1998/03/20 14:20:10 lefmann * Added a few definitions. * * Revision 1.6 1998/03/18 13:21:48 elu_krni * Added definition of PROP_Z0 macro. * * Revision 1.5 1998/03/16 08:04:16 kn * Added normal distributed random number function randnorm(). * * Revision 1.4 1997/12/03 13:34:19 kn * Added definition of ABSORB macro. * * Revision 1.3 1997/10/16 14:27:28 kn * Added debugging output used by the "display" graphical visualization * tool. * * Revision 1.2 1997/09/08 11:31:27 kn * Added mcsetstate() function. * * Revision 1.1 1997/09/08 10:39:44 kn * Initial revision * * * Copyright (C) Risoe National Laboratory, 1997-1998, All rights reserved *******************************************************************************/ #ifndef MCSTAS_R_H #define MCSTAS_R_H #include #include #include #include /* If the runtime is embedded in the simulation program, some definitions can be made static. */ #ifdef MC_EMBEDDED_RUNTIME #define mcstatic static #else #define mcstatic #endif typedef double MCNUM; typedef struct {MCNUM x, y, z;} Coords; typedef MCNUM Rotation[3][3]; struct mcinputtable_struct { char *name; MCNUM *par; }; extern struct mcinputtable_struct mcinputtable[]; extern int mcnumipar; extern char mcinstrument_name[], mcinstrument_source[]; extern int mctraceenabled, mcdefaultmain; void mcinit(void); void mcraytrace(void); void mcfinally(void); #define ABSORB do {mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz, \ mcnlt,mcnls1,mcnls2, mcnlp); mcDEBUG_ABSORB(); goto mcabsorb;} while(0) #define MC_GETPAR(comp, par) mcc ## comp ## _ ## par #define DETECTOR_OUT(p0,p1,p2) printf("Detector: %s_I=%g %s_ERR=%g\n", \ mccompcurname, p1, mccompcurname, mcestimate_error(p0,p1,p2)) #ifdef MC_TRACE_ENABLED #define DEBUG #endif #ifdef DEBUG #define mcDEBUG_INSTR() if(!mcdotrace); else printf("INSTRUMENT:\n"); #define mcDEBUG_COMPONENT(name,c,t) if(!mcdotrace); else \ printf("COMPONENT: \"%s\"\n" \ "POS: %g, %g, %g, %g, %g, %g, %g, %g, %g, %g, %g, %g\n", \ name, c.x, c.y, c.z, t[0][0], t[0][1], t[0][2], \ t[1][0], t[1][1], t[1][2], t[2][0], t[2][1], t[2][2]); #define mcDEBUG_INSTR_END() if(!mcdotrace); else printf("INSTRUMENT END:\n"); #define mcDEBUG_ENTER() if(!mcdotrace); else printf("ENTER:\n"); #define mcDEBUG_COMP(c) if(!mcdotrace); else printf("COMP: \"%s\"\n", c); #define mcDEBUG_STATE(x,y,z,vx,vy,vz,t,s1,s2,p) if(!mcdotrace); else \ printf("STATE: %g, %g, %g, %g, %g, %g, %g, %g, %g, %g\n", \ x,y,z,vx,vy,vz,t,s1,s2,p); #define mcDEBUG_LEAVE() if(!mcdotrace); else printf("LEAVE:\n"); #define mcDEBUG_ABSORB() if(!mcdotrace); else printf("ABSORB:\n"); #else #define mcDEBUG_INSTR() #define mcDEBUG_COMPONENT(name,c,t) #define mcDEBUG_INSTR_END() #define mcDEBUG_ENTER() #define mcDEBUG_COMP(c) #define mcDEBUG_STATE(x,y,z,vx,vy,vz,t,s1,s2,p) #define mcDEBUG_LEAVE() #define mcDEBUG_ABSORB() #endif #ifdef TEST #define test_printf printf #else #define test_printf while(0) printf #endif #define MIN2RAD (PI/(180*60)) #define DEG2RAD (PI/180) #define RAD2DEG (180/PI) #define AA2MS 629.719 /* Convert k[1/AA] to v[m/s] */ #define MS2AA 1.58801E-3 /* Convert v[m/s] to k[1/AA] */ #define K2V AA2MS #define V2K MS2AA #define Q2V AA2MS #define V2Q MS2AA #define SE2V 437.3949 /* Convert sqrt(E)[meV] to v[m/s] */ #define VS2E 5.227e-6 /* Convert v[m/s] to sqrt(E)[meV] */ #define HBAR 1.05459E-34 #define MNEUTRON 1.67492E-27 #ifndef PI # ifdef M_PI # define PI M_PI # else # define PI 3.14159265358979323846 # endif #endif typedef int mc_int32_t; mc_int32_t mc_random(void); void mc_srandom (unsigned int x); #ifndef USE_SYSTEM_RANDOM #ifdef RAND_MAX # undef RAND_MAX #endif #define RAND_MAX 0x7fffffff #define random mc_random #define srandom mc_srandom #endif /* !USE_SYSTEM_RANDOM */ #define rand01() ( ((double)random())/((double)RAND_MAX+1) ) #define randpm1() ( ((double)random()) / (((double)RAND_MAX+1)/2) - 1 ) #define rand0max(max) ( ((double)random()) / (((double)RAND_MAX+1)/(max)) ) #define randminmax(min,max) ( rand0max((max)-(min)) - (min) ) #define PROP_X0 \ { \ double mc_dt; \ if(mcnlvx == 0) ABSORB; \ mc_dt = -mcnlx/mcnlvx; \ if(mc_dt < 0) ABSORB; \ mcnly += mcnlvy*mc_dt; \ mcnlz += mcnlvz*mc_dt; \ mcnlt += mc_dt; \ mcnlx = 0; \ } #define PROP_Y0 \ { \ double mc_dt; \ if(mcnlvy == 0) ABSORB; \ mc_dt = -mcnly/mcnlvy; \ if(mc_dt < 0) ABSORB; \ mcnlx += mcnlvx*mc_dt; \ mcnlz += mcnlvz*mc_dt; \ mcnlt += mc_dt; \ mcnly = 0; \ } #define PROP_Z0 \ { \ double mc_dt; \ if(mcnlvz == 0) ABSORB; \ mc_dt = -mcnlz/mcnlvz; \ if(mc_dt < 0) ABSORB; \ mcnlx += mcnlvx*mc_dt; \ mcnly += mcnlvy*mc_dt; \ mcnlt += mc_dt; \ mcnlz = 0; \ } #define PROP_DT(dt) \ { \ mcnlx += mcnlvx*(dt); \ mcnly += mcnlvy*(dt); \ mcnlz += mcnlvz*(dt); \ mcnlt += (dt); \ } #define vec_prod(x, y, z, x1, y1, z1, x2, y2, z2) \ { \ double mcvp_tmpx, mcvp_tmpy, mcvp_tmpz; \ mcvp_tmpx = (y1)*(z2) - (y2)*(z1); \ mcvp_tmpy = (z1)*(x2) - (z2)*(x1); \ mcvp_tmpz = (x1)*(y2) - (x2)*(y1); \ (x) = mcvp_tmpx; (y) = mcvp_tmpy; (z) = mcvp_tmpz; \ } #define scalar_prod(x1, y1, z1, x2, y2, z2) \ ((x1)*(x2) + (y1)*(y2) + (z1)*(z2)) #define NORM(x,y,z) \ { \ double mcnm_tmp = sqrt((x)*(x) + (y)*(y) + (z)*(z)); \ if(mcnm_tmp != 0.0) \ { \ (x) /= mcnm_tmp; \ (y) /= mcnm_tmp; \ (z) /= mcnm_tmp; \ } \ } #define rotate(x, y, z, vx, vy, vz, phi, ax, ay, az) \ { \ double mcrt_tmpx = (ax), mcrt_tmpy = (ay), mcrt_tmpz = (az); \ double mcrt_vp, mcrt_vpx, mcrt_vpy, mcrt_vpz; \ double mcrt_vnx, mcrt_vny, mcrt_vnz, mcrt_vn1x, mcrt_vn1y, mcrt_vn1z; \ double mcrt_bx, mcrt_by, mcrt_bz, mcrt_v1x, mcrt_v1y, mcrt_v1z; \ double mcrt_cos, mcrt_sin; \ NORM(mcrt_tmpx, mcrt_tmpy, mcrt_tmpz); \ mcrt_vp = scalar_prod((vx), (vy), (vz), mcrt_tmpx, mcrt_tmpy, mcrt_tmpz); \ mcrt_vpx = mcrt_vp*mcrt_tmpx; \ mcrt_vpy = mcrt_vp*mcrt_tmpy; \ mcrt_vpz = mcrt_vp*mcrt_tmpz; \ mcrt_vnx = (vx) - mcrt_vpx; \ mcrt_vny = (vy) - mcrt_vpy; \ mcrt_vnz = (vz) - mcrt_vpz; \ vec_prod(mcrt_bx, mcrt_by, mcrt_bz, \ mcrt_tmpx, mcrt_tmpy, mcrt_tmpz, mcrt_vnx, mcrt_vny, mcrt_vnz); \ mcrt_cos = cos((phi)); mcrt_sin = sin((phi)); \ mcrt_vn1x = mcrt_vnx*mcrt_cos + mcrt_bx*mcrt_sin; \ mcrt_vn1y = mcrt_vny*mcrt_cos + mcrt_by*mcrt_sin; \ mcrt_vn1z = mcrt_vnz*mcrt_cos + mcrt_bz*mcrt_sin; \ (x) = mcrt_vpx + mcrt_vn1x; \ (y) = mcrt_vpy + mcrt_vn1y; \ (z) = mcrt_vpz + mcrt_vn1z; \ } Coords coords_set(MCNUM x, MCNUM y, MCNUM z); Coords coords_add(Coords a, Coords b); Coords coords_sub(Coords a, Coords b); Coords coords_neg(Coords a); void rot_set_rotation(Rotation t, double phx, double phy, double phz); void rot_mul(Rotation t1, Rotation t2, Rotation t3); void rot_copy(Rotation dest, Rotation src); void rot_transpose(Rotation src, Rotation dst); Coords rot_apply(Rotation t, Coords a); void mccoordschange(Coords a, Rotation t, double *x, double *y, double *z, double *vx, double *vy, double *vz, double *time, double *s1, double *s2); double mcestimate_error(int N, double p1, double p2); void mcreadparams(void); void mcsetstate(double x, double y, double z, double vx, double vy, double vz, double t, double s1, double s2, double p); void mcgenstate(void); double randnorm(void); int cylinder_intersect(double *t0, double *t1, double x, double y, double z, double vx, double vy, double vz, double r, double h); int sphere_intersect(double *t0, double *t1, double x, double y, double z, double vx, double vy, double vz, double r); void randvec_target_sphere(double *xo, double *yo, double *zo, double *solid_angle, double xi, double yi, double zi, double radius); void mcset_ncount(double count); int mcstas_main(int argc, char *argv[]); #endif /* MCSTAS_R_H */ /* End of file "mcstas-r.h". */ #line 335 "linup-1.c" #line 1 "mcstas-r.c" /******************************************************************************* * Runtime system for McStas. * * Project: Monte Carlo Simulation of Triple Axis Spectrometers * File name: mcstas-r.c * * Author: K.N. Aug 27, 1997 * * $Id: mcstas-r.c,v 1.15 1998/10/02 08:38:27 kn Exp kn $ * * $Log: mcstas-r.c,v $ * Revision 1.15 1998/10/02 08:38:27 kn * Added DETECTOR_OUT support. * Fixed header comment. * * Revision 1.14 1998/10/01 08:12:26 kn * Support for embedding the file in the output from McStas. * Added mcstas_main() function. * Added support for command line arguments. * * Revision 1.13 1998/09/23 13:51:35 kn * McStas now uses its own random() implementation (unless * USE_SYSTEM_RANDOM is defined). * * Revision 1.12 1998/05/19 07:54:05 kn * In randvec_target_sphere, accept radius=0 to mean that no focusing is to * be done (choose direction uniformly in full 4PI solid angle). * * Revision 1.11 1998/05/11 12:08:49 kn * Fix bug in cylinder_intersect that caused an infinite cylinder height in * some cases. * * Revision 1.10 1998/04/17 11:50:08 kn * Added sphere_intersect. * * Revision 1.9 1998/04/17 10:52:27 kn * Better names in randvec_target_sphere. * * Revision 1.8 1998/04/16 14:21:49 kn * Added randvec_target() function. * * Revision 1.7 1998/03/24 07:34:03 kn * Use rand01() in randnorm(). Fix typos. * * Revision 1.6 1998/03/20 14:19:52 lefmann * Added cylinder_intersect(). * * Revision 1.5 1998/03/16 08:03:41 kn * Added normal distributed random number function randnorm(). * * Revision 1.4 1997/10/16 14:27:05 kn * Add missing #include. Change in mcreadparams() to fit better with the * "display" visualization tool. * * Revision 1.3 1997/09/08 11:31:22 kn * Added mcsetstate() function. * * Revision 1.2 1997/09/08 11:16:43 kn * Bug fix in mccoordschange(). * * Revision 1.1 1997/09/08 10:40:03 kn * Initial revision * * * Copyright (C) Risoe National Laboratory, 1997-1998, All rights reserved *******************************************************************************/ #include #include #include #ifndef MCSTAS_R_H #include "mcstas-r.h" #endif #ifdef MC_ANCIENT_COMPATIBILITY int mctraceenabled = 0; int mcdefaultmain = 0; #endif /* Assign coordinates. */ Coords coords_set(MCNUM x, MCNUM y, MCNUM z) { Coords a; a.x = x; a.y = y; a.z = z; return a; } /* Add two coordinates. */ Coords coords_add(Coords a, Coords b) { Coords c; c.x = a.x + b.x; c.y = a.y + b.y; c.z = a.z + b.z; return c; } /* Subtract two coordinates. */ Coords coords_sub(Coords a, Coords b) { Coords c; c.x = a.x - b.x; c.y = a.y - b.y; c.z = a.z - b.z; return c; } /* Negate coordinates. */ Coords coords_neg(Coords a) { Coords b; b.x = -a.x; b.y = -a.y; b.z = -a.z; return b; } /******************************************************************************* * Get transformation for rotation first phx around x axis, then phy around y, * then phz around z. *******************************************************************************/ void rot_set_rotation(Rotation t, double phx, double phy, double phz) { double cx = cos(phx); double sx = sin(phx); double cy = cos(phy); double sy = sin(phy); double cz = cos(phz); double sz = sin(phz); t[0][0] = cy*cz; t[0][1] = sx*sy*cz + cx*sz; t[0][2] = sx*sz - cx*sy*cz; t[1][0] = -cy*sz; t[1][1] = cx*cz - sx*sy*sz; t[1][2] = sx*cz + cx*sy*sz; t[2][0] = sy; t[2][1] = -sx*cy; t[2][2] = cx*cy; } /******************************************************************************* * Matrix multiplication of transformations (this corresponds to combining * transformations). After rot_mul(T1, T2, T3), doing T3 is equal to doing * first T2, then T1. * Note that T3 must not alias (use the same array as) T1 or T2. *******************************************************************************/ void rot_mul(Rotation t1, Rotation t2, Rotation t3) { int i,j, k; for(i = 0; i < 3; i++) for(j = 0; j < 3; j++) t3[i][j] = t1[i][0]*t2[0][j] + t1[i][1]*t2[1][j] + t1[i][2]*t2[2][j]; } /******************************************************************************* * Copy a rotation transformation (needed since arrays cannot be assigned in C). *******************************************************************************/ void rot_copy(Rotation dest, Rotation src) { dest[0][0] = src[0][0]; dest[0][1] = src[0][1]; dest[0][2] = src[0][2]; dest[1][0] = src[1][0]; dest[1][1] = src[1][1]; dest[1][2] = src[1][2]; dest[2][0] = src[2][0]; dest[2][1] = src[2][1]; dest[2][2] = src[2][2]; } void rot_transpose(Rotation src, Rotation dst) { dst[0][0] = src[0][0]; dst[0][1] = src[1][0]; dst[0][2] = src[2][0]; dst[1][0] = src[0][1]; dst[1][1] = src[1][1]; dst[1][2] = src[2][1]; dst[2][0] = src[0][2]; dst[2][1] = src[1][2]; dst[2][2] = src[2][2]; } Coords rot_apply(Rotation t, Coords a) { Coords b; b.x = t[0][0]*a.x + t[0][1]*a.y + t[0][2]*a.z; b.y = t[1][0]*a.x + t[1][1]*a.y + t[1][2]*a.z; b.z = t[2][0]*a.x + t[2][1]*a.y + t[2][2]*a.z; return b; } void mccoordschange(Coords a, Rotation t, double *x, double *y, double *z, double *vx, double *vy, double *vz, double *time, double *s1, double *s2) { Coords b, c; b.x = *x; b.y = *y; b.z = *z; c = rot_apply(t, b); b = coords_add(c, a); *x = b.x; *y = b.y; *z = b.z; b.x = *vx; b.y = *vy; b.z = *vz; c = rot_apply(t, b); *vx = c.x; *vy = c.y; *vz = c.z; /* ToDo: What to do about the spin? */ } double mcestimate_error(int N, double p1, double p2) { double pmean, n1; if(N <= 1) return 0; pmean = p1 / N; n1 = N - 1; return sqrt((N/n1)*(p2 - pmean*pmean)); } void mcreadparams(void) { extern struct mcinputtable_struct mcinputtable[]; int i; for(i = 0; mcinputtable[i].name != 0; i++) { printf("Set value of instrument parameter %s:\n", mcinputtable[i].name); fflush(stdout); scanf("%lf", mcinputtable[i].par); } } void mcsetstate(double x, double y, double z, double vx, double vy, double vz, double t, double s1, double s2, double p) { extern double mcnx, mcny, mcnz, mcnvx, mcnvy, mcnvz; extern double mcnt, mcns1, mcns2, mcnp; mcnx = x; mcny = y; mcnz = z; mcnvx = vx; mcnvy = vy; mcnvz = vz; mcnt = t; mcns1 = s1; mcns2 = s2; mcnp = p; } void mcgenstate(void) { mcsetstate(0, 0, 0, 0, 0, 1, 0, 0, 0, 1); } /* McStas random number routine. */ /* * Copyright (c) 1983 Regents of the University of California. * All rights reserved. * * Redistribution and use in source and binary forms are permitted * provided that the above copyright notice and this paragraph are * duplicated in all such forms and that any documentation, * advertising materials, and other materials related to such * distribution and use acknowledge that the software was developed * by the University of California, Berkeley. The name of the * University may not be used to endorse or promote products derived * from this software without specific prior written permission. * THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR * IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED * WARRANTIES OF MERCHANTIBILITY AND FITNESS FOR A PARTICULAR PURPOSE. */ /* * This is derived from the Berkeley source: * @(#)random.c 5.5 (Berkeley) 7/6/88 * It was reworked for the GNU C Library by Roland McGrath. * Rewritten to use reentrant functions by Ulrich Drepper, 1995. */ /******************************************************************************* * Modified for McStas from glibc 2.0.7pre1 stdlib/random.c and * stdlib/random_r.c. * * This way random() is more than four times faster compared to calling * standard glibc random() on ix86 Linux, probably due to multithread support, * ELF shared library overhead, etc. It also makes McStas generated * simulations more portable (more likely to behave identically across * platforms, important for parrallel computations). *******************************************************************************/ #define TYPE_3 3 #define BREAK_3 128 #define DEG_3 31 #define SEP_3 3 static mc_int32_t randtbl[DEG_3 + 1] = { TYPE_3, -1726662223, 379960547, 1735697613, 1040273694, 1313901226, 1627687941, -179304937, -2073333483, 1780058412, -1989503057, -615974602, 344556628, 939512070, -1249116260, 1507946756, -812545463, 154635395, 1388815473, -1926676823, 525320961, -1009028674, 968117788, -123449607, 1284210865, 435012392, -2017506339, -911064859, -370259173, 1132637927, 1398500161, -205601318, }; static mc_int32_t *fptr = &randtbl[SEP_3 + 1]; static mc_int32_t *rptr = &randtbl[1]; static mc_int32_t *state = &randtbl[1]; #define rand_deg DEG_3 #define rand_sep SEP_3 static mc_int32_t *end_ptr = &randtbl[sizeof (randtbl) / sizeof (randtbl[0])]; mc_int32_t mc_random (void) { mc_int32_t result; *fptr += *rptr; /* Chucking least random bit. */ result = (*fptr >> 1) & 0x7fffffff; ++fptr; if (fptr >= end_ptr) { fptr = state; ++rptr; } else { ++rptr; if (rptr >= end_ptr) rptr = state; } return result; } void mc_srandom (unsigned int x) { /* We must make sure the seed is not 0. Take arbitrarily 1 in this case. */ state[0] = x ? x : 1; { long int i; for (i = 1; i < rand_deg; ++i) { /* This does: state[i] = (16807 * state[i - 1]) % 2147483647; but avoids overflowing 31 bits. */ long int hi = state[i - 1] / 127773; long int lo = state[i - 1] % 127773; long int test = 16807 * lo - 2836 * hi; state[i] = test + (test < 0 ? 2147483647 : 0); } fptr = &state[rand_sep]; rptr = &state[0]; for (i = 0; i < 10 * rand_deg; ++i) random (); } } /* End of McStas random number routine. */ double randnorm(void) { static double v1, v2, s; static int phase = 0; double X, u1, u2; if(phase == 0) { do { u1 = rand01(); u2 = rand01(); v1 = 2*u1 - 1; v2 = 2*u2 - 1; s = v1*v1 + v2*v2; } while(s >= 1 || s == 0); X = v1*sqrt(-2*log(s)/s); } else { X = v2*sqrt(-2*log(s)/s); } phase = 1 - phase; return X; } /* Written by: EM,NB,ABA 4.2.98 */ int cylinder_intersect(double *t0, double *t1, double x, double y, double z, double vx, double vy, double vz, double r, double h) { double v, D, t_in, t_out, y_in, y_out; v = sqrt(vx*vx+vy*vy+vz*vz); D = (2*vx*x + 2*vz*z)*(2*vx*x + 2*vz*z) - 4*(vx*vx + vz*vz)*(x*x + z*z - r*r); if (D>=0) { t_in = (-(2*vz*z + 2*vx*x) - sqrt(D))/(2*(vz*vz + vx*vx)); t_out = (-(2*vz*z + 2*vx*x) + sqrt(D))/(2*(vz*vz + vx*vx)); y_in = vy*t_in + y; y_out =vy*t_out + y; if ( (y_in > h/2 && y_out > h/2) || (y_in < -h/2 && y_out < -h/2) ) return 0; else { if (y_in > h/2) t_in = ((h/2)-y)/vy; if (y_in < -h/2) t_in = ((-h/2)-y)/vy; if (y_out > h/2) t_out = ((h/2)-y)/vy; if (y_out < -h/2) t_out = ((-h/2)-y)/vy; } *t0 = t_in; *t1 = t_out; return 1; } else { *t0 = *t1 = 0; return 0; } } /* Calculate intersection between line and sphere. */ int sphere_intersect(double *t0, double *t1, double x, double y, double z, double vx, double vy, double vz, double r) { double A, B, C, D, v; v = sqrt(vx*vx + vy*vy + vz*vz); A = v*v; B = 2*(x*vx + y*vy + z*vz); C = x*x + y*y + z*z - r*r; D = B*B - 4*A*C; if(D < 0) return 0; D = sqrt(D); *t0 = (-B - D) / (2*A); *t1 = (-B + D) / (2*A); return 1; } /* Choose random direction towards target at (x,y,z) with given radius. */ /* If radius is zero, choose random direction in full 4PI, no target. */ /* ToDo: It should be possible to optimize this to avoid computing angles. */ void randvec_target_sphere(double *xo, double *yo, double *zo, double *solid_angle, double xi, double yi, double zi, double radius) { double l, theta0, phi, theta, nx, ny, nz, xt, yt, zt; if(radius == 0.0) { /* No target, choose uniformly a direction in full 4PI solid angle. */ theta = acos (1 - rand0max(2)); phi = rand0max(2 * PI); if(solid_angle) *solid_angle = 4*PI; nx = 1; ny = 0; nz = 0; xi = 0; yi = 1; zi = 0; } else { l = sqrt(xi*xi + yi*yi + zi*zi); /* Distance to target. */ theta0 = asin(radius / l); /* Target size as seen from origin */ if(solid_angle) { /* Compute solid angle of target as seen from origin. */ *solid_angle = 2*PI*(1 - cos(theta0)); } /* Now choose point uniformly on sphere surface within angle theta0 */ theta = acos (1 - rand0max(1 - cos(theta0))); phi = rand0max(2 * PI); /* Now, to obtain the desired vector rotate (x,y,z) angle theta around a perpendicular axis (nx,ny,nz) and then angle phi around (x,y,z). */ if(xi == 0 && yi == 0) { nx = 1; ny = 0; nz = 0; } else { nx = -yi; ny = xi; nz = 0; } } rotate(xt, yt, zt, xi, yi, zi, theta, nx, ny, nz); rotate(*xo, *yo, *zo, xt, yt, zt, phi, xi, yi, zi); } /* Number of neutron histories to simulate. */ static double mcncount = 1e6; void mcset_ncount(double count) { mcncount = count; } mcstatic int mcdotrace = 0; static void mcsetn_arg(char *arg) { mcset_ncount(strtod(arg, NULL)); } static void mcsetseed(char *arg) { srandom(atol(arg)); } static void mchelp(char *pgmname) { int i; fprintf(stderr, "Usage: %s [options] [parm=value ...]\n", pgmname); fprintf(stderr, "Options are:\n" " -s seed --seed=seed Set random seed\n" " -n count --ncount=count Set number of neutrons to simulate.\n" " -t --trace Enable trace of neutron through instrument\n" " -h --help Show help message\n" " -i --info Detailed instrument information\n" "Instrument parameters are:\n"); for(i = 0; i < mcnumipar; i++) fprintf(stderr, " %s\n", mcinputtable[i].name); } static void mcshowhelp(char *pgmname) { mchelp(pgmname); exit(0); } static void mcusage(char *pgmname) { fprintf(stderr, "Error: incorrect commad line arguments\n"); mchelp(pgmname); exit(1); } static void mcinfo(void) { int i; printf("Name: %s\n", mcinstrument_name); printf("Parameters:"); for(i = 0; i < mcnumipar; i++) printf(" %s", mcinputtable[i].name); printf("\n"); printf("Instrument-source: %s\n", mcinstrument_source); printf("Trace-enabled: %s\n", mctraceenabled ? "yes" : "no"); printf("Default-main: %s\n", mcdefaultmain ? "yes" : "no"); printf("Embedded-runtime: %s\n", #ifdef MC_EMBEDDED_RUNTIME "yes" #else "no" #endif ); exit(0); } static void mcenabletrace(void) { if(mctraceenabled) mcdotrace = 1; else { fprintf(stderr, "Error: trace not enabled.\n" "Please re-run the McStas compiler with the --trace option\n"); exit(1); } } void mcparseoptions(int argc, char *argv[]) { int i, j, pos; char *p; int paramset = 0, *paramsetarray; paramsetarray = malloc(mcnumipar*sizeof(*paramsetarray)); if(paramsetarray == NULL) { fprintf(stderr, "Error: insufficient memory\n"); exit(1); } for(j = 0; j < mcnumipar; j++) paramsetarray[j] = 0; for(i = 1; i < argc; i++) { if(!strcmp("-s", argv[i]) && (i + 1) < argc) mcsetseed(argv[++i]); else if(!strncmp("-s", argv[i], 2)) mcsetseed(&argv[i][2]); else if(!strcmp("--seed", argv[i]) && (i + 1) < argc) mcsetseed(argv[++i]); else if(!strncmp("--seed=", argv[i], 7)) mcsetseed(&argv[i][7]); else if(!strcmp("-n", argv[i]) && (i + 1) < argc) mcsetn_arg(argv[++i]); else if(!strncmp("-n", argv[i], 2)) mcsetn_arg(&argv[i][2]); else if(!strcmp("--ncount", argv[i]) && (i + 1) < argc) mcsetn_arg(argv[++i]); else if(!strncmp("--ncount=", argv[i], 9)) mcsetn_arg(&argv[i][9]); else if(!strcmp("-h", argv[i])) mcshowhelp(argv[0]); else if(!strcmp("--help", argv[i])) mcshowhelp(argv[0]); else if(!strcmp("-i", argv[i])) mcinfo(); else if(!strcmp("--info", argv[i])) mcinfo(); else if(!strcmp("-t", argv[i])) mcenabletrace(); else if(!strcmp("--trace", argv[i])) mcenabletrace(); else if(argv[i][0] != '-' && (p = strchr(argv[i], '=')) != NULL) { *p++ = '\0'; for(j = 0; j < mcnumipar; j++) if(!strcmp(mcinputtable[j].name, argv[i])) { *(mcinputtable[j].par) = strtod(p, NULL); paramsetarray[j] = 1; paramset = 1; break; } if(j == mcnumipar) { /* Unrecognized parameter name */ fprintf(stderr, "Error: unrecognized parameter %s\n", argv[i]); exit(1); } } else mcusage(argv[0]); } if(!paramset) mcreadparams(); /* Prompt for parameters if not specified. */ else { for(j = 0; j < mcnumipar; j++) if(!paramsetarray[j]) { fprintf(stderr, "Error: Instrument parameter %s left unset\n", mcinputtable[j].name); exit(1); } } } /* McStas main() function. */ int mcstas_main(int argc, char *argv[]) { int run_num = 0; srandom(time(NULL)); mcparseoptions(argc, argv); mcinit(); while(run_num < mcncount) { mcsetstate(0, 0, 0, 0, 0, 1, 0, 0, 0, 1); mcraytrace(); run_num++; } mcfinally(); return 0; } /* End of file "mcstas-r.c". */ #line 1084 "linup-1.c" #ifdef MC_TRACE_ENABLED int mctraceenabled = 1; #else int mctraceenabled = 0; #endif int mcdefaultmain = 1; char mcinstrument_name[] = "TAS1"; char mcinstrument_source[] = "linup-1.instr"; int main(int argc, char *argv[]){return mcstas_main(argc, argv);} void mcinit(void); void mcraytrace(void); void mcfinally(void); void mcdisplay(void); /* Instrument parameters. */ MCNUM mcipPHM; MCNUM mcipTTM; MCNUM mcipC1; #define mcNUMIPAR 3 int mcnumipar = 3; struct mcinputtable_struct mcinputtable[mcNUMIPAR+1] = { "PHM", &mcipPHM, "TTM", &mcipTTM, "C1", &mcipC1, NULL, NULL }; /* User declarations from instrument definition. */ #define PHM mcipPHM #define TTM mcipTTM #define C1 mcipC1 #line 5 "linup-1.instr" /* Mosaicity used on monochromator and analysator */ double tas1_mono_mosaic = 45; /* Measurements indicate its really 45' */ /* Q vector for bragg scattering with monochromator and analysator */ double tas1_mono_q = 3.354; /* Fake 2nd order scattering for 20meV */ double tas1_mono_r0 = 0.6; double mpos0, mpos1, mpos2, mpos3, mpos4, mpos5, mpos6, mpos7; double mrot0, mrot1, mrot2, mrot3, mrot4, mrot5, mrot6, mrot7; #line 1126 "linup-1.c" #undef C1 #undef TTM #undef PHM /* Declarations of component definition and setting parameters. */ /* Definition parameters for component 'source'. */ #define mccsource_radius 0.06 #define mccsource_dist 3.288 #define mccsource_xw 0.042 #define mccsource_yh 0.082 #define mccsource_E0 20 #define mccsource_dE 0.82 /* Definition parameters for component 'slit1'. */ #define mccslit1_xmin -0.02 #define mccslit1_xmax 0.065 #define mccslit1_ymin -0.075 #define mccslit1_ymax 0.075 /* Definition parameters for component 'slit2'. */ #define mccslit2_xmin -0.02 #define mccslit2_xmax 0.02 #define mccslit2_ymin -0.04 #define mccslit2_ymax 0.04 /* Definition parameters for component 'slit3'. */ #define mccslit3_xmin -0.021 #define mccslit3_xmax 0.021 #define mccslit3_ymin -0.041 #define mccslit3_ymax 0.041 /* Definition parameters for component 'm0'. */ #define mccm0_zmin -0.0375 #define mccm0_zmax 0.0375 #define mccm0_ymin -0.006 #define mccm0_ymax 0.006 #define mccm0_mosaich tas1_mono_mosaic #define mccm0_mosaicv tas1_mono_mosaic #define mccm0_r0 tas1_mono_r0 #define mccm0_Q tas1_mono_q /* Definition parameters for component 'm1'. */ #define mccm1_zmin -0.0375 #define mccm1_zmax 0.0375 #define mccm1_ymin -0.006 #define mccm1_ymax 0.006 #define mccm1_mosaich tas1_mono_mosaic #define mccm1_mosaicv tas1_mono_mosaic #define mccm1_r0 tas1_mono_r0 #define mccm1_Q tas1_mono_q /* Definition parameters for component 'm2'. */ #define mccm2_zmin -0.0375 #define mccm2_zmax 0.0375 #define mccm2_ymin -0.006 #define mccm2_ymax 0.006 #define mccm2_mosaich tas1_mono_mosaic #define mccm2_mosaicv tas1_mono_mosaic #define mccm2_r0 tas1_mono_r0 #define mccm2_Q tas1_mono_q /* Definition parameters for component 'm3'. */ #define mccm3_zmin -0.0375 #define mccm3_zmax 0.0375 #define mccm3_ymin -0.006 #define mccm3_ymax 0.006 #define mccm3_mosaich tas1_mono_mosaic #define mccm3_mosaicv tas1_mono_mosaic #define mccm3_r0 tas1_mono_r0 #define mccm3_Q tas1_mono_q /* Definition parameters for component 'm4'. */ #define mccm4_zmin -0.0375 #define mccm4_zmax 0.0375 #define mccm4_ymin -0.006 #define mccm4_ymax 0.006 #define mccm4_mosaich tas1_mono_mosaic #define mccm4_mosaicv tas1_mono_mosaic #define mccm4_r0 tas1_mono_r0 #define mccm4_Q tas1_mono_q /* Definition parameters for component 'm5'. */ #define mccm5_zmin -0.0375 #define mccm5_zmax 0.0375 #define mccm5_ymin -0.006 #define mccm5_ymax 0.006 #define mccm5_mosaich tas1_mono_mosaic #define mccm5_mosaicv tas1_mono_mosaic #define mccm5_r0 tas1_mono_r0 #define mccm5_Q tas1_mono_q /* Definition parameters for component 'm6'. */ #define mccm6_zmin -0.0375 #define mccm6_zmax 0.0375 #define mccm6_ymin -0.006 #define mccm6_ymax 0.006 #define mccm6_mosaich tas1_mono_mosaic #define mccm6_mosaicv tas1_mono_mosaic #define mccm6_r0 tas1_mono_r0 #define mccm6_Q tas1_mono_q /* Definition parameters for component 'm7'. */ #define mccm7_zmin -0.0375 #define mccm7_zmax 0.0375 #define mccm7_ymin -0.006 #define mccm7_ymax 0.006 #define mccm7_mosaich tas1_mono_mosaic #define mccm7_mosaicv tas1_mono_mosaic #define mccm7_r0 tas1_mono_r0 #define mccm7_Q tas1_mono_q /* Definition parameters for component 'slitMS1'. */ #define mccslitMS1_xmin -0.0105 #define mccslitMS1_xmax 0.0105 #define mccslitMS1_ymin -0.035 #define mccslitMS1_ymax 0.035 /* Definition parameters for component 'slitMS2'. */ #define mccslitMS2_xmin -0.0105 #define mccslitMS2_xmax 0.0105 #define mccslitMS2_ymin -0.035 #define mccslitMS2_ymax 0.035 /* Definition parameters for component 'c1'. */ #define mccc1_xmin -0.02 #define mccc1_xmax 0.02 #define mccc1_ymin -0.0375 #define mccc1_ymax 0.0375 #define mccc1_len 0.25 #define mccc1_divergence mcipC1 /* Definition parameters for component 'slitMS3'. */ #define mccslitMS3_radius 0.025 /* Definition parameters for component 'slitMS4'. */ #define mccslitMS4_radius 0.025 /* Definition parameters for component 'slitMS5'. */ #define mccslitMS5_radius 0.0275 /* Definition parameters for component 'emon1'. */ #define mccemon1_xmin -0.01 #define mccemon1_xmax 0.01 #define mccemon1_ymin -0.1 #define mccemon1_ymax 0.1 #define mccemon1_Emin 19.25 #define mccemon1_Emax 20.75 #define mccemon1_nchan 35 #define mccemon1_filename "sim/linup_1_1.emon" /* Definition parameters for component 'slitS'. */ #define mccslitS_xmin -0.00525 #define mccslitS_xmax 0.00525 #define mccslitS_ymin -0.02025 #define mccslitS_ymax 0.02025 /* Definition parameters for component 'sng'. */ #define mccsng_xmin -0.025 #define mccsng_xmax 0.025 #define mccsng_ymin -0.0375 #define mccsng_ymax 0.0375 /* User component declarations. */ /* User declarations for component 'source'. */ #define mccompcurname "source" #define radius mccsource_radius #define dist mccsource_dist #define xw mccsource_xw #define yh mccsource_yh #define E0 mccsource_E0 #define dE mccsource_dE #define hdiv mccsource_hdiv #define vdiv mccsource_vdiv #define p_in mccsource_p_in #line 36 "/data/lnslib/lib/mcstas/Source_flat.comp" double hdiv,vdiv; double p_in; #line 1306 "linup-1.c" #undef p_in #undef vdiv #undef hdiv #undef dE #undef E0 #undef yh #undef xw #undef dist #undef radius #undef mccompcurname /* User declarations for component 'm0'. */ #define mccompcurname "m0" #define zmin mccm0_zmin #define zmax mccm0_zmax #define ymin mccm0_ymin #define ymax mccm0_ymax #define mosaich mccm0_mosaich #define mosaicv mccm0_mosaicv #define r0 mccm0_r0 #define Q mccm0_Q #line 37 "/data/lnslib/lib/mcstas/Monochromator.comp" #define DIV_CUTOFF 2 /* ~ 10^-5 cutoff. */ #line 1330 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname /* User declarations for component 'm1'. */ #define mccompcurname "m1" #define zmin mccm1_zmin #define zmax mccm1_zmax #define ymin mccm1_ymin #define ymax mccm1_ymax #define mosaich mccm1_mosaich #define mosaicv mccm1_mosaicv #define r0 mccm1_r0 #define Q mccm1_Q #line 37 "/data/lnslib/lib/mcstas/Monochromator.comp" #define DIV_CUTOFF 2 /* ~ 10^-5 cutoff. */ #line 1353 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname /* User declarations for component 'm2'. */ #define mccompcurname "m2" #define zmin mccm2_zmin #define zmax mccm2_zmax #define ymin mccm2_ymin #define ymax mccm2_ymax #define mosaich mccm2_mosaich #define mosaicv mccm2_mosaicv #define r0 mccm2_r0 #define Q mccm2_Q #line 37 "/data/lnslib/lib/mcstas/Monochromator.comp" #define DIV_CUTOFF 2 /* ~ 10^-5 cutoff. */ #line 1376 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname /* User declarations for component 'm3'. */ #define mccompcurname "m3" #define zmin mccm3_zmin #define zmax mccm3_zmax #define ymin mccm3_ymin #define ymax mccm3_ymax #define mosaich mccm3_mosaich #define mosaicv mccm3_mosaicv #define r0 mccm3_r0 #define Q mccm3_Q #line 37 "/data/lnslib/lib/mcstas/Monochromator.comp" #define DIV_CUTOFF 2 /* ~ 10^-5 cutoff. */ #line 1399 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname /* User declarations for component 'm4'. */ #define mccompcurname "m4" #define zmin mccm4_zmin #define zmax mccm4_zmax #define ymin mccm4_ymin #define ymax mccm4_ymax #define mosaich mccm4_mosaich #define mosaicv mccm4_mosaicv #define r0 mccm4_r0 #define Q mccm4_Q #line 37 "/data/lnslib/lib/mcstas/Monochromator.comp" #define DIV_CUTOFF 2 /* ~ 10^-5 cutoff. */ #line 1422 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname /* User declarations for component 'm5'. */ #define mccompcurname "m5" #define zmin mccm5_zmin #define zmax mccm5_zmax #define ymin mccm5_ymin #define ymax mccm5_ymax #define mosaich mccm5_mosaich #define mosaicv mccm5_mosaicv #define r0 mccm5_r0 #define Q mccm5_Q #line 37 "/data/lnslib/lib/mcstas/Monochromator.comp" #define DIV_CUTOFF 2 /* ~ 10^-5 cutoff. */ #line 1445 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname /* User declarations for component 'm6'. */ #define mccompcurname "m6" #define zmin mccm6_zmin #define zmax mccm6_zmax #define ymin mccm6_ymin #define ymax mccm6_ymax #define mosaich mccm6_mosaich #define mosaicv mccm6_mosaicv #define r0 mccm6_r0 #define Q mccm6_Q #line 37 "/data/lnslib/lib/mcstas/Monochromator.comp" #define DIV_CUTOFF 2 /* ~ 10^-5 cutoff. */ #line 1468 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname /* User declarations for component 'm7'. */ #define mccompcurname "m7" #define zmin mccm7_zmin #define zmax mccm7_zmax #define ymin mccm7_ymin #define ymax mccm7_ymax #define mosaich mccm7_mosaich #define mosaicv mccm7_mosaicv #define r0 mccm7_r0 #define Q mccm7_Q #line 37 "/data/lnslib/lib/mcstas/Monochromator.comp" #define DIV_CUTOFF 2 /* ~ 10^-5 cutoff. */ #line 1491 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname /* User declarations for component 'c1'. */ #define mccompcurname "c1" #define xmin mccc1_xmin #define xmax mccc1_xmax #define ymin mccc1_ymin #define ymax mccc1_ymax #define len mccc1_len #define divergence mccc1_divergence #define slope mccc1_slope #line 36 "/data/lnslib/lib/mcstas/Soller.comp" double slope; #line 1513 "linup-1.c" #undef slope #undef divergence #undef len #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname /* User declarations for component 'emon1'. */ #define mccompcurname "emon1" #define xmin mccemon1_xmin #define xmax mccemon1_xmax #define ymin mccemon1_ymin #define ymax mccemon1_ymax #define Emin mccemon1_Emin #define Emax mccemon1_Emax #define nchan mccemon1_nchan #define filename mccemon1_filename #define E_N mccemon1_E_N #define E_p mccemon1_E_p #define E_p2 mccemon1_E_p2 #line 40 "/data/lnslib/lib/mcstas/E_monitor.comp" int E_N[nchan]; double E_p[nchan], E_p2[nchan]; #line 1539 "linup-1.c" #undef E_p2 #undef E_p #undef E_N #undef filename #undef nchan #undef Emax #undef Emin #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname /* User declarations for component 'sng'. */ #define mccompcurname "sng" #define xmin mccsng_xmin #define xmax mccsng_xmax #define ymin mccsng_ymin #define ymax mccsng_ymax #define Nsum mccsng_Nsum #define psum mccsng_psum #define p2sum mccsng_p2sum #line 36 "/data/lnslib/lib/mcstas/Monitor.comp" int Nsum; double psum, p2sum; #line 1565 "linup-1.c" #undef p2sum #undef psum #undef Nsum #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname Coords mcposaa1, mcposra1; Rotation mcrotaa1, mcrotra1; Coords mcposasource, mcposrsource; Rotation mcrotasource, mcrotrsource; Coords mcposaslit1, mcposrslit1; Rotation mcrotaslit1, mcrotrslit1; Coords mcposaslit2, mcposrslit2; Rotation mcrotaslit2, mcrotrslit2; Coords mcposaslit3, mcposrslit3; Rotation mcrotaslit3, mcrotrslit3; Coords mcposafocus_mono, mcposrfocus_mono; Rotation mcrotafocus_mono, mcrotrfocus_mono; Coords mcposam0, mcposrm0; Rotation mcrotam0, mcrotrm0; Coords mcposam1, mcposrm1; Rotation mcrotam1, mcrotrm1; Coords mcposam2, mcposrm2; Rotation mcrotam2, mcrotrm2; Coords mcposam3, mcposrm3; Rotation mcrotam3, mcrotrm3; Coords mcposam4, mcposrm4; Rotation mcrotam4, mcrotrm4; Coords mcposam5, mcposrm5; Rotation mcrotam5, mcrotrm5; Coords mcposam6, mcposrm6; Rotation mcrotam6, mcrotrm6; Coords mcposam7, mcposrm7; Rotation mcrotam7, mcrotrm7; Coords mcposaa2, mcposra2; Rotation mcrotaa2, mcrotra2; Coords mcposaslitMS1, mcposrslitMS1; Rotation mcrotaslitMS1, mcrotrslitMS1; Coords mcposaslitMS2, mcposrslitMS2; Rotation mcrotaslitMS2, mcrotrslitMS2; Coords mcposac1, mcposrc1; Rotation mcrotac1, mcrotrc1; Coords mcposaslitMS3, mcposrslitMS3; Rotation mcrotaslitMS3, mcrotrslitMS3; Coords mcposaslitMS4, mcposrslitMS4; Rotation mcrotaslitMS4, mcrotrslitMS4; Coords mcposaslitMS5, mcposrslitMS5; Rotation mcrotaslitMS5, mcrotrslitMS5; Coords mcposaemon1, mcposremon1; Rotation mcrotaemon1, mcrotremon1; Coords mcposaa3, mcposra3; Rotation mcrotaa3, mcrotra3; Coords mcposaslitS, mcposrslitS; Rotation mcrotaslitS, mcrotrslitS; Coords mcposasng, mcposrsng; Rotation mcrotasng, mcrotrsng; MCNUM mcnx, mcny, mcnz, mcnvx, mcnvy, mcnvz, mcnt, mcnsx, mcnsy, mcnsz, mcnp; void mcinit(void) { #define PHM mcipPHM #define TTM mcipTTM #define C1 mcipC1 #line 16 "linup-1.instr" { double d = 0.0125; /* 12.5 mm between slab centers. */ double phi = 0.5443; /* Rotation between adjacent slabs. */ mpos0 = -3.5*d; mrot0 = -3.5*phi; mpos1 = -2.5*d; mrot1 = -2.5*phi; mpos2 = -1.5*d; mrot2 = -1.5*phi; mpos3 = -0.5*d; mrot3 = -0.5*phi; mpos4 = 0.5*d; mrot4 = 0.5*phi; mpos5 = 1.5*d; mrot5 = 1.5*phi; mpos6 = 2.5*d; mrot6 = 2.5*phi; mpos7 = 3.5*d; mrot7 = 3.5*phi; } #line 1645 "linup-1.c" #undef C1 #undef TTM #undef PHM /* Component initializations. */ /* Initializations for component a1. */ /* Initializations for component source. */ #define mccompcurname "source" #define radius mccsource_radius #define dist mccsource_dist #define xw mccsource_xw #define yh mccsource_yh #define E0 mccsource_E0 #define dE mccsource_dE #define hdiv mccsource_hdiv #define vdiv mccsource_vdiv #define p_in mccsource_p_in #line 40 "/data/lnslib/lib/mcstas/Source_flat.comp" { hdiv = atan(xw/(2.0*dist)); vdiv = atan(yh/(2.0*dist)); p_in = (4*hdiv*vdiv)/(4*PI); /* Small angle approx. */ } #line 1671 "linup-1.c" #undef p_in #undef vdiv #undef hdiv #undef dE #undef E0 #undef yh #undef xw #undef dist #undef radius #undef mccompcurname /* Initializations for component slit1. */ /* Initializations for component slit2. */ /* Initializations for component slit3. */ /* Initializations for component focus_mono. */ /* Initializations for component m0. */ /* Initializations for component m1. */ /* Initializations for component m2. */ /* Initializations for component m3. */ /* Initializations for component m4. */ /* Initializations for component m5. */ /* Initializations for component m6. */ /* Initializations for component m7. */ /* Initializations for component a2. */ /* Initializations for component slitMS1. */ /* Initializations for component slitMS2. */ /* Initializations for component c1. */ #define mccompcurname "c1" #define xmin mccc1_xmin #define xmax mccc1_xmax #define ymin mccc1_ymin #define ymax mccc1_ymax #define len mccc1_len #define divergence mccc1_divergence #define slope mccc1_slope #line 39 "/data/lnslib/lib/mcstas/Soller.comp" { slope = tan(MIN2RAD*divergence); } #line 1742 "linup-1.c" #undef slope #undef divergence #undef len #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname /* Initializations for component slitMS3. */ /* Initializations for component slitMS4. */ /* Initializations for component slitMS5. */ /* Initializations for component emon1. */ #define mccompcurname "emon1" #define xmin mccemon1_xmin #define xmax mccemon1_xmax #define ymin mccemon1_ymin #define ymax mccemon1_ymax #define Emin mccemon1_Emin #define Emax mccemon1_Emax #define nchan mccemon1_nchan #define filename mccemon1_filename #define E_N mccemon1_E_N #define E_p mccemon1_E_p #define E_p2 mccemon1_E_p2 #line 44 "/data/lnslib/lib/mcstas/E_monitor.comp" { int i; for (i=0; ixmax || yymax) ABSORB; } #line 2269 "linup-1.c" #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component slit2. */ mcDEBUG_COMP("slit2") mccoordschange(mcposrslit2, mcrotrslit2, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "slit2" #define xmin mccslit2_xmin #define xmax mccslit2_xmax #define ymin mccslit2_ymin #define ymax mccslit2_ymax #line 28 "/data/lnslib/lib/mcstas/Slit.comp" { PROP_Z0; if (xxmax || yymax) ABSORB; } #line 2315 "linup-1.c" #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component slit3. */ mcDEBUG_COMP("slit3") mccoordschange(mcposrslit3, mcrotrslit3, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "slit3" #define xmin mccslit3_xmin #define xmax mccslit3_xmax #define ymin mccslit3_ymin #define ymax mccslit3_ymax #line 28 "/data/lnslib/lib/mcstas/Slit.comp" { PROP_Z0; if (xxmax || yymax) ABSORB; } #line 2361 "linup-1.c" #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component focus_mono. */ mcDEBUG_COMP("focus_mono") mccoordschange(mcposrfocus_mono, mcrotrfocus_mono, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define mccompcurname "focus_mono" #undef mccompcurname mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component m0. */ mcDEBUG_COMP("m0") mccoordschange(mcposrm0, mcrotrm0, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "m0" #define zmin mccm0_zmin #define zmax mccm0_zmax #define ymin mccm0_ymin #define ymax mccm0_ymax #define mosaich mccm0_mosaich #define mosaicv mccm0_mosaicv #define r0 mccm0_r0 #define Q mccm0_Q #line 41 "/data/lnslib/lib/mcstas/Monochromator.comp" { double dphi,tmp1,tmp2,tmp3,vratio,phi,theta0,theta,v,cs,sn; double old_x = x, old_y = y, old_z = z, old_t = t; double dt; if(vx != 0.0 && (dt = -x/vx) >= 0.0) { y += vy*dt; z += vz*dt; t += dt; x = 0.0; if (z>zmin && zymin && y DIV_CUTOFF) { x = old_x; y = old_y; z = old_z; t = old_t; } else { p *= r0*exp(-tmp3*tmp3*4*log(2)); /* Use mosaics */ tmp1 = 2*theta; cs = cos(tmp1); sn = sin(tmp1); tmp2 = cs*vx - sn*vz; vy = vy; vz = cs*vz + sn*vx; vx = tmp2; /* Second: scatering out of plane. Approximation is that Debye-Scherrer cone is a plane */ phi = atan2(vy,vz); /* out-of plane angle */ dphi = (MIN2RAD*mosaicv)/(2*sqrt(2*log(2)))*randnorm(); /* MC choice: */ /* Vertical angle of the crystallite */ vy = vz*tan(phi+2*dphi*sin(theta)); vratio = v/sqrt(vx*vx+vy*vy+vz*vz); vz = vz*vratio; vy = vy*vratio; /* Renormalize v */ vx = vx*vratio; } } else { x = old_x; y = old_y; z = old_z; t = old_t; } } } #line 2470 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component m1. */ mcDEBUG_COMP("m1") mccoordschange(mcposrm1, mcrotrm1, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "m1" #define zmin mccm1_zmin #define zmax mccm1_zmax #define ymin mccm1_ymin #define ymax mccm1_ymax #define mosaich mccm1_mosaich #define mosaicv mccm1_mosaicv #define r0 mccm1_r0 #define Q mccm1_Q #line 41 "/data/lnslib/lib/mcstas/Monochromator.comp" { double dphi,tmp1,tmp2,tmp3,vratio,phi,theta0,theta,v,cs,sn; double old_x = x, old_y = y, old_z = z, old_t = t; double dt; if(vx != 0.0 && (dt = -x/vx) >= 0.0) { y += vy*dt; z += vz*dt; t += dt; x = 0.0; if (z>zmin && zymin && y DIV_CUTOFF) { x = old_x; y = old_y; z = old_z; t = old_t; } else { p *= r0*exp(-tmp3*tmp3*4*log(2)); /* Use mosaics */ tmp1 = 2*theta; cs = cos(tmp1); sn = sin(tmp1); tmp2 = cs*vx - sn*vz; vy = vy; vz = cs*vz + sn*vx; vx = tmp2; /* Second: scatering out of plane. Approximation is that Debye-Scherrer cone is a plane */ phi = atan2(vy,vz); /* out-of plane angle */ dphi = (MIN2RAD*mosaicv)/(2*sqrt(2*log(2)))*randnorm(); /* MC choice: */ /* Vertical angle of the crystallite */ vy = vz*tan(phi+2*dphi*sin(theta)); vratio = v/sqrt(vx*vx+vy*vy+vz*vz); vz = vz*vratio; vy = vy*vratio; /* Renormalize v */ vx = vx*vratio; } } else { x = old_x; y = old_y; z = old_z; t = old_t; } } } #line 2572 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component m2. */ mcDEBUG_COMP("m2") mccoordschange(mcposrm2, mcrotrm2, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "m2" #define zmin mccm2_zmin #define zmax mccm2_zmax #define ymin mccm2_ymin #define ymax mccm2_ymax #define mosaich mccm2_mosaich #define mosaicv mccm2_mosaicv #define r0 mccm2_r0 #define Q mccm2_Q #line 41 "/data/lnslib/lib/mcstas/Monochromator.comp" { double dphi,tmp1,tmp2,tmp3,vratio,phi,theta0,theta,v,cs,sn; double old_x = x, old_y = y, old_z = z, old_t = t; double dt; if(vx != 0.0 && (dt = -x/vx) >= 0.0) { y += vy*dt; z += vz*dt; t += dt; x = 0.0; if (z>zmin && zymin && y DIV_CUTOFF) { x = old_x; y = old_y; z = old_z; t = old_t; } else { p *= r0*exp(-tmp3*tmp3*4*log(2)); /* Use mosaics */ tmp1 = 2*theta; cs = cos(tmp1); sn = sin(tmp1); tmp2 = cs*vx - sn*vz; vy = vy; vz = cs*vz + sn*vx; vx = tmp2; /* Second: scatering out of plane. Approximation is that Debye-Scherrer cone is a plane */ phi = atan2(vy,vz); /* out-of plane angle */ dphi = (MIN2RAD*mosaicv)/(2*sqrt(2*log(2)))*randnorm(); /* MC choice: */ /* Vertical angle of the crystallite */ vy = vz*tan(phi+2*dphi*sin(theta)); vratio = v/sqrt(vx*vx+vy*vy+vz*vz); vz = vz*vratio; vy = vy*vratio; /* Renormalize v */ vx = vx*vratio; } } else { x = old_x; y = old_y; z = old_z; t = old_t; } } } #line 2674 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component m3. */ mcDEBUG_COMP("m3") mccoordschange(mcposrm3, mcrotrm3, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "m3" #define zmin mccm3_zmin #define zmax mccm3_zmax #define ymin mccm3_ymin #define ymax mccm3_ymax #define mosaich mccm3_mosaich #define mosaicv mccm3_mosaicv #define r0 mccm3_r0 #define Q mccm3_Q #line 41 "/data/lnslib/lib/mcstas/Monochromator.comp" { double dphi,tmp1,tmp2,tmp3,vratio,phi,theta0,theta,v,cs,sn; double old_x = x, old_y = y, old_z = z, old_t = t; double dt; if(vx != 0.0 && (dt = -x/vx) >= 0.0) { y += vy*dt; z += vz*dt; t += dt; x = 0.0; if (z>zmin && zymin && y DIV_CUTOFF) { x = old_x; y = old_y; z = old_z; t = old_t; } else { p *= r0*exp(-tmp3*tmp3*4*log(2)); /* Use mosaics */ tmp1 = 2*theta; cs = cos(tmp1); sn = sin(tmp1); tmp2 = cs*vx - sn*vz; vy = vy; vz = cs*vz + sn*vx; vx = tmp2; /* Second: scatering out of plane. Approximation is that Debye-Scherrer cone is a plane */ phi = atan2(vy,vz); /* out-of plane angle */ dphi = (MIN2RAD*mosaicv)/(2*sqrt(2*log(2)))*randnorm(); /* MC choice: */ /* Vertical angle of the crystallite */ vy = vz*tan(phi+2*dphi*sin(theta)); vratio = v/sqrt(vx*vx+vy*vy+vz*vz); vz = vz*vratio; vy = vy*vratio; /* Renormalize v */ vx = vx*vratio; } } else { x = old_x; y = old_y; z = old_z; t = old_t; } } } #line 2776 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component m4. */ mcDEBUG_COMP("m4") mccoordschange(mcposrm4, mcrotrm4, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "m4" #define zmin mccm4_zmin #define zmax mccm4_zmax #define ymin mccm4_ymin #define ymax mccm4_ymax #define mosaich mccm4_mosaich #define mosaicv mccm4_mosaicv #define r0 mccm4_r0 #define Q mccm4_Q #line 41 "/data/lnslib/lib/mcstas/Monochromator.comp" { double dphi,tmp1,tmp2,tmp3,vratio,phi,theta0,theta,v,cs,sn; double old_x = x, old_y = y, old_z = z, old_t = t; double dt; if(vx != 0.0 && (dt = -x/vx) >= 0.0) { y += vy*dt; z += vz*dt; t += dt; x = 0.0; if (z>zmin && zymin && y DIV_CUTOFF) { x = old_x; y = old_y; z = old_z; t = old_t; } else { p *= r0*exp(-tmp3*tmp3*4*log(2)); /* Use mosaics */ tmp1 = 2*theta; cs = cos(tmp1); sn = sin(tmp1); tmp2 = cs*vx - sn*vz; vy = vy; vz = cs*vz + sn*vx; vx = tmp2; /* Second: scatering out of plane. Approximation is that Debye-Scherrer cone is a plane */ phi = atan2(vy,vz); /* out-of plane angle */ dphi = (MIN2RAD*mosaicv)/(2*sqrt(2*log(2)))*randnorm(); /* MC choice: */ /* Vertical angle of the crystallite */ vy = vz*tan(phi+2*dphi*sin(theta)); vratio = v/sqrt(vx*vx+vy*vy+vz*vz); vz = vz*vratio; vy = vy*vratio; /* Renormalize v */ vx = vx*vratio; } } else { x = old_x; y = old_y; z = old_z; t = old_t; } } } #line 2878 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component m5. */ mcDEBUG_COMP("m5") mccoordschange(mcposrm5, mcrotrm5, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "m5" #define zmin mccm5_zmin #define zmax mccm5_zmax #define ymin mccm5_ymin #define ymax mccm5_ymax #define mosaich mccm5_mosaich #define mosaicv mccm5_mosaicv #define r0 mccm5_r0 #define Q mccm5_Q #line 41 "/data/lnslib/lib/mcstas/Monochromator.comp" { double dphi,tmp1,tmp2,tmp3,vratio,phi,theta0,theta,v,cs,sn; double old_x = x, old_y = y, old_z = z, old_t = t; double dt; if(vx != 0.0 && (dt = -x/vx) >= 0.0) { y += vy*dt; z += vz*dt; t += dt; x = 0.0; if (z>zmin && zymin && y DIV_CUTOFF) { x = old_x; y = old_y; z = old_z; t = old_t; } else { p *= r0*exp(-tmp3*tmp3*4*log(2)); /* Use mosaics */ tmp1 = 2*theta; cs = cos(tmp1); sn = sin(tmp1); tmp2 = cs*vx - sn*vz; vy = vy; vz = cs*vz + sn*vx; vx = tmp2; /* Second: scatering out of plane. Approximation is that Debye-Scherrer cone is a plane */ phi = atan2(vy,vz); /* out-of plane angle */ dphi = (MIN2RAD*mosaicv)/(2*sqrt(2*log(2)))*randnorm(); /* MC choice: */ /* Vertical angle of the crystallite */ vy = vz*tan(phi+2*dphi*sin(theta)); vratio = v/sqrt(vx*vx+vy*vy+vz*vz); vz = vz*vratio; vy = vy*vratio; /* Renormalize v */ vx = vx*vratio; } } else { x = old_x; y = old_y; z = old_z; t = old_t; } } } #line 2980 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component m6. */ mcDEBUG_COMP("m6") mccoordschange(mcposrm6, mcrotrm6, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "m6" #define zmin mccm6_zmin #define zmax mccm6_zmax #define ymin mccm6_ymin #define ymax mccm6_ymax #define mosaich mccm6_mosaich #define mosaicv mccm6_mosaicv #define r0 mccm6_r0 #define Q mccm6_Q #line 41 "/data/lnslib/lib/mcstas/Monochromator.comp" { double dphi,tmp1,tmp2,tmp3,vratio,phi,theta0,theta,v,cs,sn; double old_x = x, old_y = y, old_z = z, old_t = t; double dt; if(vx != 0.0 && (dt = -x/vx) >= 0.0) { y += vy*dt; z += vz*dt; t += dt; x = 0.0; if (z>zmin && zymin && y DIV_CUTOFF) { x = old_x; y = old_y; z = old_z; t = old_t; } else { p *= r0*exp(-tmp3*tmp3*4*log(2)); /* Use mosaics */ tmp1 = 2*theta; cs = cos(tmp1); sn = sin(tmp1); tmp2 = cs*vx - sn*vz; vy = vy; vz = cs*vz + sn*vx; vx = tmp2; /* Second: scatering out of plane. Approximation is that Debye-Scherrer cone is a plane */ phi = atan2(vy,vz); /* out-of plane angle */ dphi = (MIN2RAD*mosaicv)/(2*sqrt(2*log(2)))*randnorm(); /* MC choice: */ /* Vertical angle of the crystallite */ vy = vz*tan(phi+2*dphi*sin(theta)); vratio = v/sqrt(vx*vx+vy*vy+vz*vz); vz = vz*vratio; vy = vy*vratio; /* Renormalize v */ vx = vx*vratio; } } else { x = old_x; y = old_y; z = old_z; t = old_t; } } } #line 3082 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component m7. */ mcDEBUG_COMP("m7") mccoordschange(mcposrm7, mcrotrm7, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "m7" #define zmin mccm7_zmin #define zmax mccm7_zmax #define ymin mccm7_ymin #define ymax mccm7_ymax #define mosaich mccm7_mosaich #define mosaicv mccm7_mosaicv #define r0 mccm7_r0 #define Q mccm7_Q #line 41 "/data/lnslib/lib/mcstas/Monochromator.comp" { double dphi,tmp1,tmp2,tmp3,vratio,phi,theta0,theta,v,cs,sn; double old_x = x, old_y = y, old_z = z, old_t = t; double dt; if(vx != 0.0 && (dt = -x/vx) >= 0.0) { y += vy*dt; z += vz*dt; t += dt; x = 0.0; if (z>zmin && zymin && y DIV_CUTOFF) { x = old_x; y = old_y; z = old_z; t = old_t; } else { p *= r0*exp(-tmp3*tmp3*4*log(2)); /* Use mosaics */ tmp1 = 2*theta; cs = cos(tmp1); sn = sin(tmp1); tmp2 = cs*vx - sn*vz; vy = vy; vz = cs*vz + sn*vx; vx = tmp2; /* Second: scatering out of plane. Approximation is that Debye-Scherrer cone is a plane */ phi = atan2(vy,vz); /* out-of plane angle */ dphi = (MIN2RAD*mosaicv)/(2*sqrt(2*log(2)))*randnorm(); /* MC choice: */ /* Vertical angle of the crystallite */ vy = vz*tan(phi+2*dphi*sin(theta)); vratio = v/sqrt(vx*vx+vy*vy+vz*vz); vz = vz*vratio; vy = vy*vratio; /* Renormalize v */ vx = vx*vratio; } } else { x = old_x; y = old_y; z = old_z; t = old_t; } } } #line 3184 "linup-1.c" #undef Q #undef r0 #undef mosaicv #undef mosaich #undef ymax #undef ymin #undef zmax #undef zmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component a2. */ mcDEBUG_COMP("a2") mccoordschange(mcposra2, mcrotra2, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define mccompcurname "a2" #undef mccompcurname mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component slitMS1. */ mcDEBUG_COMP("slitMS1") mccoordschange(mcposrslitMS1, mcrotrslitMS1, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "slitMS1" #define xmin mccslitMS1_xmin #define xmax mccslitMS1_xmax #define ymin mccslitMS1_ymin #define ymax mccslitMS1_ymax #line 28 "/data/lnslib/lib/mcstas/Slit.comp" { PROP_Z0; if (xxmax || yymax) ABSORB; } #line 3245 "linup-1.c" #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component slitMS2. */ mcDEBUG_COMP("slitMS2") mccoordschange(mcposrslitMS2, mcrotrslitMS2, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "slitMS2" #define xmin mccslitMS2_xmin #define xmax mccslitMS2_xmax #define ymin mccslitMS2_ymin #define ymax mccslitMS2_ymax #line 28 "/data/lnslib/lib/mcstas/Slit.comp" { PROP_Z0; if (xxmax || yymax) ABSORB; } #line 3291 "linup-1.c" #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component c1. */ mcDEBUG_COMP("c1") mccoordschange(mcposrc1, mcrotrc1, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "c1" #define xmin mccc1_xmin #define xmax mccc1_xmax #define ymin mccc1_ymin #define ymax mccc1_ymax #define len mccc1_len #define divergence mccc1_divergence #define slope mccc1_slope #line 43 "/data/lnslib/lib/mcstas/Soller.comp" { double phi, dt; PROP_Z0; if (xxmax || yymax) ABSORB; dt = len/vz; PROP_DT(dt); if (xxmax || yymax) ABSORB; if(slope > 0.0) { phi = fabs(vx/vz); if (phi > slope) ABSORB; else p *= 1.0 - phi/slope; } } #line 3355 "linup-1.c" #undef slope #undef divergence #undef len #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component slitMS3. */ mcDEBUG_COMP("slitMS3") mccoordschange(mcposrslitMS3, mcrotrslitMS3, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "slitMS3" #define radius mccslitMS3_radius #line 22 "/data/lnslib/lib/mcstas/Circular_slit.comp" { PROP_Z0; if(x*x + y*y > radius*radius) ABSORB; } #line 3401 "linup-1.c" #undef radius #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component slitMS4. */ mcDEBUG_COMP("slitMS4") mccoordschange(mcposrslitMS4, mcrotrslitMS4, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "slitMS4" #define radius mccslitMS4_radius #line 22 "/data/lnslib/lib/mcstas/Circular_slit.comp" { PROP_Z0; if(x*x + y*y > radius*radius) ABSORB; } #line 3441 "linup-1.c" #undef radius #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component slitMS5. */ mcDEBUG_COMP("slitMS5") mccoordschange(mcposrslitMS5, mcrotrslitMS5, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "slitMS5" #define radius mccslitMS5_radius #line 22 "/data/lnslib/lib/mcstas/Circular_slit.comp" { PROP_Z0; if(x*x + y*y > radius*radius) ABSORB; } #line 3481 "linup-1.c" #undef radius #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component emon1. */ mcDEBUG_COMP("emon1") mccoordschange(mcposremon1, mcrotremon1, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "emon1" #define xmin mccemon1_xmin #define xmax mccemon1_xmax #define ymin mccemon1_ymin #define ymax mccemon1_ymax #define Emin mccemon1_Emin #define Emax mccemon1_Emax #define nchan mccemon1_nchan #define filename mccemon1_filename #define E_N mccemon1_E_N #define E_p mccemon1_E_p #define E_p2 mccemon1_E_p2 #line 55 "/data/lnslib/lib/mcstas/E_monitor.comp" { int i; double E; PROP_Z0; if (x>xmin && xymin && y= 0 && i < nchan) { E_N[i]++; E_p[i] += p; E_p2[i] += p*p; } } } #line 3544 "linup-1.c" #undef E_p2 #undef E_p #undef E_N #undef filename #undef nchan #undef Emax #undef Emin #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component a3. */ mcDEBUG_COMP("a3") mccoordschange(mcposra3, mcrotra3, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define mccompcurname "a3" #undef mccompcurname mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component slitS. */ mcDEBUG_COMP("slitS") mccoordschange(mcposrslitS, mcrotrslitS, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "slitS" #define xmin mccslitS_xmin #define xmax mccslitS_xmax #define ymin mccslitS_ymin #define ymax mccslitS_ymax #line 28 "/data/lnslib/lib/mcstas/Slit.comp" { PROP_Z0; if (xxmax || yymax) ABSORB; } #line 3608 "linup-1.c" #undef ymax #undef ymin #undef xmax #undef xmin #undef mccompcurname #undef p #undef s2 #undef s1 #undef t #undef vz #undef vy #undef vx #undef z #undef y #undef x mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) /* Component sng. */ mcDEBUG_COMP("sng") mccoordschange(mcposrsng, mcrotrsng, &mcnlx, &mcnly, &mcnlz, &mcnlvx, &mcnlvy, &mcnlvz, &mcnlt, &mcnlsx, &mcnlsy); mcDEBUG_STATE(mcnlx, mcnly, mcnlz, mcnlvx, mcnlvy, mcnlvz,mcnlt,mcnlsx,mcnlsy, mcnlp) #define x mcnlx #define y mcnly #define z mcnlz #define vx mcnlvx #define vy mcnlvy #define vz mcnlvz #define t mcnlt #define s1 mcnlsx #define s2 mcnlsy #define p mcnlp #define mccompcurname "sng" #define xmin mccsng_xmin #define xmax mccsng_xmax #define ymin mccsng_ymin #define ymax mccsng_ymax #define Nsum mccsng_Nsum #define psum mccsng_psum #define p2sum mccsng_p2sum #line 46 "/data/lnslib/lib/mcstas/Monitor.comp" { PROP_Z0; if (x>xmin && xymin && y (Bertrand.Roessli@psi.ch) Message-ID: > Dear Kristian, > > please find below the linup-1.c file which > has got undefined variables. Ok, from this I was able to spot the cause of the problems with undefined symbols. The problem is the mixing of two different versions of McStas. The file has been generated with version 1.1 of McStas, but using files from version 1.0. This means that the 'mcstas' command that was run is from version 1.1, but when the program searched for the files it needed (particularly mcstas-r.c and mcstas-r.h), it found the wrong ones from version 1.0. This will not work. I wonder how this has happened? Perhaps the most likely cause is that when you installed version 1.0, the environment variable MCSTAS was set to point to the location of the 1.0 files (you can check this with the command echo $MCSTAS and see if it gives the name of a directory). Then McStas v1.1 was installed in another location, without updating the environment variable. Usually, the environment variable is not necessary since McStas remembers the location of its files when it is installed. Setting the variable to the wrong value would cause the problems you experienced. I hope this clarifies things for you? Otherwise, perhaps you can write me the exact command that was used to produce the file 'linup-1.c' along with the output from the "echo $MCSTAS" command above? - Kristian. From Georg_Artus at Physik.TU-Muenchen.DE Fri Jun 25 11:10:42 1999 From: Georg_Artus at Physik.TU-Muenchen.DE (Dr. Georg Artus) Date: Fri, 25 Jun 1999 11:10:42 +0200 Subject: Soller_trans.comp References: <01JCJPKXZQ7C94NTAD@risoe.dk> Message-ID: <37734792.D83D6BB9@ph.tum.de> Hello Kristian, I just modified the component Soller to Soller_trans. The new component has an additional parameter transmission, which is a general multiplicator for the transmission probability. The true transmission of a Soller is not 1 also with a non divergent beam. For calculations of absolute fluxes this is important. Because I am not experience with C (I just checked Soller.comp by mathematical undterstanding) I don?t want to send the new component to the mailing-list directly. I also didn?t add a check whether the given parameter ?transmission? lies between 0 and 1. Our tests indicate that the component works correctly. May be you can have a look at the file. Georg -- ********************************************* Dr. Georg Artus Technische Universitaet Muenchen Neue Forschungs-Neutronenquelle Garching D-85747 Garching Tel: +49 (0)89/289-14675 Fax: +49 (0)89/289-14666 or 12112 E-mail: gartus at ph.tum.de -------------- next part -------------- /******************************************************************************* * * McStas, version 1.0, released October 26, 1998 * Maintained by Kristian Nielsen and Kim Lefmann, * Risoe National Laboratory, Roskilde, Denmark * * Component: Soller_trans * * Written by: KN, August 1998 * Modified by: GA, June 1999 * * Soller collimator with rectangular opening and specified length. The * transmission function is an average and does not utilize knowledge of the * actual neutron trajectory. * A zero divergence disables collimation (then the component works as a double slit). * Added is an additional parameter ?transmission? to take into account the * transmission probability. * * INPUT PARAMETERS: * * xmin: (m) Lower x bound on slits * xmax: (m) Upper x bound on slits * ymin: (m) Lower y bound on slits * ymax: (m) Upper y bound on slits * len: (m) Distance between slits * divergence: (minutes of arc) Divergence angle (calculated as atan(d/len), * where d is the blade spacing) * transmission: (1) 0<=transmission<=1) * *******************************************************************************/ DEFINE COMPONENT Soller_trans DEFINITION PARAMETERS (xmin, xmax, ymin, ymax, len, divergence, transmission) SETTING PARAMETERS () OUTPUT PARAMETERS (slope) STATE PARAMETERS (x,y,z,vx,vy,vz,t,s1,s2,p) DECLARE %{ double slope; %} INITIALIZE %{ slope = tan(MIN2RAD*divergence); %} TRACE %{ double phi, dt; PROP_Z0; if (xxmax || yymax) ABSORB; dt = len/vz; PROP_DT(dt); if (xxmax || yymax) ABSORB; if(slope > 0.0) { phi = fabs(vx/vz); if (phi > slope) ABSORB; else p *= transmission*(1.0 - phi/slope); } %} MCDISPLAY %{ double x; int i; magnify("xy"); for(x = xmin, i = 0; i <= 3; i++, x += (xmax - xmin)/3.0) multiline(5, x, (double)ymin, 0.0, x, (double)ymax, 0.0, x, (double)ymax, (double)len, x, (double)ymin, (double)len, x, (double)ymin, 0.0); line(xmin, ymin, 0, xmax, ymin, 0); line(xmin, ymax, 0, xmax, ymax, 0); line(xmin, ymin, len, xmax, ymin, len); line(xmin, ymax, len, xmax, ymax, len); %} END From Stephan_Roth at Physik.TU-Muenchen.DE Wed Jun 23 09:04:20 1999 From: Stephan_Roth at Physik.TU-Muenchen.DE (Stephan Roth) Date: Wed, 23 Jun 1999 09:04:20 +0200 (MET DST) Subject: Vanadium Sample Message-ID: Dear Sir, I am using McStas for the construction of a disc chopper time-of-flight spectrometer of IN5-type. Especially I want to simulate angular resolved time-of-flight spectra at the detector bank. That is I want to see the influence of the finite diameter of cylinder samples on the spectra at several scattering angles. So I want to use a focusing option like in V_sample.comp, but the focusing should be on a rectangular detector. So I first made a setup of IN5-type with a Ni-guide, a V-sample and a realistic secondary flight path of 4m, in order to take a look the angular spectrum at a scattering angle of 90degree (without tof-detector). But a detector placed at the focusing point counts no neutrons. Hence I would like to ask you, if you could help me finding the mistake in the setup (I have incorporated the in5.instr file in this mail). Thank you in advance for your help. Best regards Stephan Roth _/_/_/_/ _/_/ _/_/_/_/ Stephan V. Roth _/ _/ _/ _/ _/ Technische Universitaet Muenchen _/ _/ _/ _/ Physik-Department Lehrstuhl E13 _/_/_/ _/ _/_/ James-Franck-Str. 1 _/ _/ _/ 85748 Garching _/ _/ _/ _/ Tel.: +49 (0)89-289-12883 _/_/_/_/ _/ _/_/_/_/ Fax : +49 (0)89-289-12473 -------------- next part -------------- DEFINE INSTRUMENT detektest() TRACE COMPONENT arm = Arm() AT (0,0,0) ABSOLUTE COMPONENT source = Source_div(width=.02,height=.05,hdiv=.702,vdiv=.702,E0=2.2722,dE=.001) AT (0,0,0) RELATIVE arm COMPONENT guide = Guide(w1=.02,h1=0.05,w2=.02,h2=0.05,l=53,R0=0.995, Qc=.0214,alpha=5.566,m=1,W=.0033) AT (0,0,0.001) RELATIVE source COMPONENT detector1=Monitor(xmin=-.1,xmax=.1,ymin=-.2,ymax=.2) AT (0,0,0.2) RELATIVE guide COMPONENT sample1=V_sample(radius_i=.0143,radius_o=.0148,h=.05,pack=1,focus_r=.1, target_x=4,target_y=0,target_z=0) AT (0,0,.205) RELATIVE guide /* COMPONENT slit_det=Circular_slit(radius=.1) AT (4.000,0,0) RELATIVE sample1 */ COMPONENT detector1f=Monitor(xmin=-.2,xmax=.2,ymin=-.2,ymax=.2) AT (4,0,0) RELATIVE sample1 END