[3] | 1 | /*
|
---|
| 2 | ---- Hector the simulator ----
|
---|
| 3 | A fast simulator of particles through generic beamlines.
|
---|
| 4 | J. de Favereau, X. Rouby ~~~ hector_devel@cp3.phys.ucl.ac.be
|
---|
| 5 |
|
---|
| 6 | http://www.fynu.ucl.ac.be/hector.html
|
---|
| 7 |
|
---|
| 8 | Centre de Physique des Particules et de Phénoménologie (CP3)
|
---|
| 9 | Université Catholique de Louvain (UCL)
|
---|
| 10 | */
|
---|
| 11 |
|
---|
| 12 | /// \file H_RectEllipticAperture.cc
|
---|
| 13 | /// \brief Defines the Rect-Elliptic aperture of beamline elements.
|
---|
| 14 |
|
---|
| 15 | // C++ #includes
|
---|
| 16 | #include <iostream>
|
---|
| 17 |
|
---|
| 18 | // C #includes
|
---|
| 19 | #include <cmath> // needed for fabs
|
---|
| 20 |
|
---|
| 21 | // ROOT #includes
|
---|
| 22 | #include "TPolyLine.h"
|
---|
| 23 |
|
---|
| 24 | // local #includes
|
---|
| 25 | #include "H_RectEllipticAperture.h"
|
---|
| 26 | using namespace std;
|
---|
| 27 |
|
---|
| 28 | H_RectEllipticAperture::H_RectEllipticAperture(const float l, const float h, const float L, const float H, const float posx, const float posy) :H_Aperture(RECTELLIPSE,((l==0)?L:l),((h==0)?H:h),L,H,posx,posy) {
|
---|
| 29 | /// @param l, h, L, H are the geometrical parameters of the rect-ellipse shape
|
---|
| 30 | /// @param posx, posy defines the (x,y) of the center of the shape
|
---|
| 31 | type= RECTELLIPSE;
|
---|
| 32 | }
|
---|
| 33 |
|
---|
| 34 | H_RectEllipticAperture* H_RectEllipticAperture::clone() const {
|
---|
| 35 | return new H_RectEllipticAperture(x1,x2,x3,x4,fx,fy);
|
---|
| 36 | }
|
---|
| 37 |
|
---|
| 38 |
|
---|
| 39 | TPolyLine * rectellipse(const float a_e = 2, const float b_e = 1, const float a_r = 1, const float b_r = 2, const float center_x = 0, const float center_y =0) {
|
---|
| 40 | const int n = 20; // number of points per segment
|
---|
| 41 | const int N = 4*n; // there are 4 segments
|
---|
| 42 | float x[N+1], y[N+1];
|
---|
| 43 |
|
---|
| 44 | if(a_e>a_r) {
|
---|
| 45 | // a rectellipse has 4 segments
|
---|
| 46 | // 1) upper one
|
---|
| 47 | for (int i=0; i<n; i++) {
|
---|
| 48 | x[i] = -a_r + i*(2*a_r)/(float)n;
|
---|
| 49 | y[i] = b_e * sqrt(1-pow(x[i]/a_e,2));
|
---|
| 50 | }
|
---|
| 51 |
|
---|
| 52 | // 2) right vertical segment
|
---|
| 53 | // upper right corner
|
---|
| 54 | const float y2 = b_e * sqrt(1-pow(a_r/a_e,2));
|
---|
| 55 | // lower right corner
|
---|
| 56 | const float y3 = -b_e * sqrt(1-pow(a_r/a_e,2));
|
---|
| 57 | for (int i=n; i<2*n; i++) {
|
---|
| 58 | x[i] = a_r;
|
---|
| 59 | y[i] = y2 - (i-n)*(2*y2)/(float)n;
|
---|
| 60 | }
|
---|
| 61 |
|
---|
| 62 | // 3) lower side
|
---|
| 63 | for (int i=2*n; i<3*n; i++) {
|
---|
| 64 | x[i] = a_r - (i-2*n)*(2*a_r)/(float)n;
|
---|
| 65 | y[i] = -b_e * sqrt(1-pow(x[i]/a_e,2));
|
---|
| 66 | }
|
---|
| 67 |
|
---|
| 68 | // 4) left vertical segment
|
---|
| 69 | // lower left corner
|
---|
| 70 | const float y4 = y3;
|
---|
| 71 | for (int i=3*n; i<4*n; i++) {
|
---|
| 72 | x[i] = -a_r;
|
---|
| 73 | y[i] = y4 + (i-3*n)*(2*y2)/(float)n;
|
---|
| 74 | }
|
---|
| 75 | } else {
|
---|
| 76 | // 1) upper one : flat
|
---|
| 77 | const float x1 = -a_e * sqrt(1-pow(b_r/b_e,2));
|
---|
| 78 | const float x2 = -x1;
|
---|
| 79 | for (int i=0; i<n; i++) {
|
---|
| 80 | y[i] = b_r;
|
---|
| 81 | x[i] = x1 + i * (x2-x1)/(float)n;
|
---|
| 82 | }
|
---|
| 83 |
|
---|
| 84 | // 2) right curved border
|
---|
| 85 | for (int i=n; i<2*n; i++) {
|
---|
| 86 | y[i] = b_r - (i-n) * (2*b_r)/(float)n;
|
---|
| 87 | x[i] = a_e * sqrt(1-pow(y[i]/b_e,2));
|
---|
| 88 | }
|
---|
| 89 |
|
---|
| 90 | // 3) lower side : flat
|
---|
| 91 | for (int i=2*n; i<3*n; i++) {
|
---|
| 92 | y[i] = -b_r;
|
---|
| 93 | x[i] = x2 - (i-2*n) * (2*x2)/(float)n;
|
---|
| 94 | }
|
---|
| 95 |
|
---|
| 96 | // 4) left curved border
|
---|
| 97 | for (int i=3*n; i<4*n; i++) {
|
---|
| 98 | y[i] = -b_r + (i-3*n) * (2*b_r)/(float)n;
|
---|
| 99 | x[i] = -a_e * sqrt(1-pow(y[i]/b_e,2));
|
---|
| 100 | }
|
---|
| 101 | }
|
---|
| 102 |
|
---|
| 103 | // closing the polyline
|
---|
| 104 | x[N] = x[0];
|
---|
| 105 | y[N] = y[0];
|
---|
| 106 |
|
---|
| 107 | // shifting the center
|
---|
| 108 | for (int i=0; i<N+1; i++) {
|
---|
| 109 | x[i] += center_x;
|
---|
| 110 | y[i] += center_y;
|
---|
| 111 | }
|
---|
| 112 |
|
---|
| 113 | return new TPolyLine(N+1,x,y);
|
---|
| 114 | }
|
---|
| 115 |
|
---|
| 116 |
|
---|
| 117 | void H_RectEllipticAperture::draw(const float scale) const {
|
---|
| 118 | TPolyLine * re = rectellipse(x3*scale, x4*scale, x1*scale, x2*scale, fx*scale, fy*scale);
|
---|
| 119 | re->SetLineColor(39);
|
---|
| 120 | re->SetLineWidth(2);
|
---|
| 121 | re->Draw("l");
|
---|
| 122 | return;
|
---|
| 123 | }
|
---|
| 124 |
|
---|
| 125 | bool H_RectEllipticAperture::isInside(const float x, const float y) const {
|
---|
| 126 | return((fabs(fx-x)<x1)&&(fabs(fy-y)<x2)&&(((x-fx)/x3)*((x-fx)/x3)+((y-fy)/x4)*((y-fy)/x4)<1));
|
---|
| 127 | }
|
---|
| 128 |
|
---|
| 129 | void H_RectEllipticAperture::printProperties() const {
|
---|
| 130 | cout << "Aperture shape:" << getTypeString() << ", parameters " <<x1<<", "<<x2<<", "<<x3<<", "<<x4<< endl;
|
---|
| 131 | cout << " \t Center : "<<fx<<", "<<fy<<endl;
|
---|
| 132 | return;
|
---|
| 133 | }
|
---|