diff --git a/cpp/src/core/sroptcryst.h b/cpp/src/core/sroptcryst.h index 19a6607c..13064348 100644 --- a/cpp/src/core/sroptcryst.h +++ b/cpp/src/core/sroptcryst.h @@ -10,6 +10,8 @@ * * @author J.Sutter, A.Suvorov, O.Chubar * @version 1.0 + * @author A. Rodriguez-Fernandez, I. Petrov + * @version 2.0 #Laue ***************************************************************************/ #ifndef __SROPTCRYST_H @@ -52,12 +54,14 @@ class srTOptCryst : public srTGenOptElem { //TVector3d m_nv; // horizontal, vertical and longitudinal coordinates of outward normal to crystal surface in the frame of incident beam TVector3d m_tv; // horizontal, vertical and longitudinal coordinates of central tangential vector [m] in the frame of incident beam TVector3d m_sv; // horizontal, vertical and longitudinal coordinates of central saggital vector [m] in the frame of incident beam + TVector3d v_nv; //* horizontal, vertical and longitudinal coordinates of outward normal to crystal surface in the frame of incident beam */ //ARF06062019 Introduce Vector3d& v_nv, - char m_uc; // crystal use case: 1- Bragg Reflection, 2- Bragg Transmission - // Laue cases added: 3- Laue, 4- Laue Transmission ARF08052019 - - // TRANSFORMATION MATRIX - double m_RXtLab[3][3]; // RXtLab: 3x3 orthogonal matrix that converts components of a 3x1 vector in crystal coordinates to components in lab coordinates. + int m_uc; // crystal use case: 1- Bragg Reflection, 2- Bragg Transmission (Laue cases to be added) + // case: 3- Laue, 4- Laue Transmission ARF08052019 change char to int + int m_ug; //ARF2305219 This variable define Laue (2) or Bragg (1) + + // TRANSFORMATIN MATRIX + double m_RXtLab[3][3]; // RXtLab: 3x3 orthogonal matrix that converts components of a 3x1 vector in crystal coordinates to components in lab coordinates. double m_RLabXt[3][3]; // RLabXt: transpose of RXtLab: converts components of a 3x1 vector in lab coordinates to components in crystal coordinates. double m_RXtRef[3][3]; // RXtRef: 3x3 orthogonal matrix that converts components of a 3x1 vector in crystal coordinates to components in diffracted beam coordinates. @@ -73,7 +77,9 @@ class srTOptCryst : public srTGenOptElem { double m_sg0X[3]; double m_cos2t; //cos(2.*thBrd); - + double thBrd; //ARF23052019 to calculate the gamma0 and gammaH + double alphrd; //ARF23052019 to calculate the gamma0 and gammaH + srTOptCrystMeshTrf *m_pMeshTrf; double m_eStartAux, m_eStepAux, m_ne; //bool m_xStartIsConstInSlicesE, m_zStartIsConstInSlicesE; @@ -82,7 +88,7 @@ class srTOptCryst : public srTGenOptElem { srTOptCryst(const SRWLOptCryst& srwlCr) { - m_dA = srwlCr.dSp; //crystal reflecting planes d-spacing (units: [A] ?) + m_dA = srwlCr.dSp; //crystal reflecting planes d-spacing (units: [A] ) //m_psi0r = srwlCr.psi0r; //m_psi0i = srwlCr.psi0i; @@ -98,53 +104,54 @@ class srTOptCryst : public srTGenOptElem { //m_kMilND = srwlCr.h2; //m_lMilND = srwlCr.h3; //OC180314: the Miller indices are removed after discussion with A. Suvorov (because these are only required for the m_dA, and it is used as input parameter) - //m_tc = srwlCr.tc; - m_thicum = srwlCr.tc*1e+06; //assuming srwlCr.tc in [m] - - //ARF22053019 Geometry - m_uc = srwlCr.uc; //OC05092016 - if((m_uc < 1) || (m_uc > 4)) throw IMPROPER_OPTICAL_COMPONENT_CRYSTAL_USE_CASE; //ARF 08052019 Change m_uc>2 to 4 for adding Laue cases - + + //Thickness units [microns] ? + m_thicum = srwlCr.tc * 1e+06; //assuming srwlCr.tc in [m] goes to [um] + + //Geometry + m_uc = (int) srwlCr.uc; //OC05092016 + if((m_uc < 1) || (m_uc > 4)) throw IMPROPER_OPTICAL_COMPONENT_CRYSTAL_USE_CASE; //ARF 08052019 Change m_uc>2 to 4 for adding Laue cases. + // Input #7: Whether to calculate the transmitted beam as well as the diffracted + // beam. itrans = 0 for NO, itrans = 1 for YES. + m_ug = 1; //Bragg case ARF22053019 + if ((m_uc == 3) || (m_uc == 4)) m_ug = 2; //Laue case ARF22053019 + m_itrans = 0; //OC: make it input variable - if((m_uc == 2) || (m_uc == 4)) m_itrans = 1; //ARF05082019 add m_uc == 4 - - double alphrd = srwlCr.angAs; //asymmetry angle [rad] - + if ((m_uc == 2) || (m_uc == 4)) m_itrans = 1; //OC05092016 //ARF05082019 add m_uc == 4 + + + alphrd = srwlCr.angAs; //asymmetry angle [rad] + // From Input #5: reciprocal lattice vector coordinates - //m_HXAi[0] = 0.; - //m_HXAi[1] = cos(alphrd) / m_dA; - //m_HXAi[2] = -sin(alphrd) / m_dA; - //OC02102019, following ARF - if((m_uc == 1) || (m_uc == 2)) - {//Bragg - m_HXAi[0] = 0.; - m_HXAi[1] = cos(alphrd)/m_dA; - m_HXAi[2] = -sin(alphrd)/m_dA; - } - else if((m_uc == 3) || (m_uc == 4)) - {//Laue - //IP14062019 change recipr. lat. vector for Laue + m_HXAi[0] = 0.; + m_HXAi[1] = cos(alphrd) / m_dA; + m_HXAi[2] = -sin(alphrd) / m_dA; + if (m_ug == 2) + { m_HXAi[0] = 0.; - m_HXAi[1] = -sin(alphrd)/m_dA; + m_HXAi[1] = sin(alphrd)/m_dA; m_HXAi[2] = -cos(alphrd)/m_dA; } - if(srwlCr.nvz == 0) throw IMPROPER_OPTICAL_COMPONENT_ORIENT; - TVector3d v_nv(srwlCr.nvx, srwlCr.nvy, srwlCr.nvz); /* horizontal, vertical and longitudinal coordinates of outward normal to crystal surface in the frame of incident beam */ - v_nv.Normalize(); - double nv[] = { v_nv.x, v_nv.y, v_nv.z }; - + v_nv.x = srwlCr.nvx; v_nv.y = srwlCr.nvy; v_nv.z = srwlCr.nvz; + v_nv.Normalize(); + double nv[] = { v_nv.x, v_nv.y, v_nv.z }; //ARF change from double nv[] = { v_nv.x, v_nv.y, v_nv.z }; + v_nv.x = nv[0];v_nv.y = nv[1];v_nv.z = nv[2]; //ARF06062019 Introduce Vector3d& v_nv, + if((srwlCr.tvx == 0) && (srwlCr.tvy == 0)) throw IMPROPER_OPTICAL_COMPONENT_ORIENT; //TVector3d v_tv(srwlCr.tvx, srwlCr.tvy, 0); /* horizontal, vertical and vertical coordinates of central tangential vector [m] in the frame of incident beam */ //v_tv.z = (-v_nv.x*v_tv.x - v_nv.y*v_tv.y) / v_nv.z; //v_tv.Normalize(); //double tv[] = { v_tv.x, v_tv.y, v_tv.z }; - m_tv.x = srwlCr.tvx; m_tv.y = srwlCr.tvy; /* horizontal and vertical coordinates of central tangential vector [m] in the frame of incident beam */ + + m_tv.x = srwlCr.tvx; m_tv.y = srwlCr.tvy; /* horizontal and vertical coordinates of central tangential vector [m] in the frame of incident beam */ m_tv.z = (-v_nv.x*m_tv.x - v_nv.y*m_tv.y) / v_nv.z; + m_tv.Normalize(); + double tv[] = { m_tv.x, m_tv.y, m_tv.z }; //sv[0] = nv[1] * tv[2] - nv[2] * tv[1]; @@ -153,14 +160,6 @@ class srTOptCryst : public srTGenOptElem { double sv[] = { nv[1] * tv[2] - nv[2] * tv[1], nv[2] * tv[0] - nv[0] * tv[2], nv[0] * tv[1] - nv[1] * tv[0] }; m_sv.x = sv[0]; m_sv.y = sv[1]; m_sv.z = sv[2]; - //m_uc = srwlCr.uc; //OC02102019 (moved up, as suggested by ARF) //OC05092016 - //if((m_uc < 1) || (m_uc > 2)) throw IMPROPER_OPTICAL_COMPONENT_MIRROR_USE_CASE; - - // Input #7: Whether to calculate the transmitted beam as well as the diffracted - // beam. itrans = 0 for NO, itrans = 1 for YES. - //OC02102019 (moved up, as suggested by ARF) - //m_itrans = 0; //OC: make it input variable - //if(m_uc == 2) m_itrans = 1; //OC05092016 // RXtLab: 3x3 orthogonal matrix that converts components of a 3x1 vector // in crystal coordinates to components in lab coordinates. @@ -174,30 +173,15 @@ class srTOptCryst : public srTGenOptElem { // coordinates to components in crystal coordinates. for(int i = 0; i < 3; i++) for(int j = 0; j < 3; j++) m_RLabXt[i][j] = m_RXtLab[j][i]; - // Dynamical Bragg angle (angle between incident optical axis and lattice planes) double uBrd = m_RLabXt[0][2] * m_HXAi[0] + m_RLabXt[1][2] * m_HXAi[1] + m_RLabXt[2][2] * m_HXAi[2]; double uHX = sqrt(m_HXAi[0] * m_HXAi[0] + m_HXAi[1] * m_HXAi[1] + m_HXAi[2] * m_HXAi[2]); double uRL = sqrt(m_RLabXt[0][2] * m_RLabXt[0][2] + m_RLabXt[1][2] * m_RLabXt[1][2] + m_RLabXt[2][2] * m_RLabXt[2][2]); - const double pi = 4.*atan(1.); - double thBrd = acos(uBrd / uHX / uRL) - pi / 2.; + + thBrd = acos(uBrd / uHX / uRL) - pi / 2.; //printf("Bragg angle = %f\n",thBrd*180./pi); m_cos2t = cos(2.*thBrd); - //ARF03062019 - if((m_uc == 1) || (m_uc == 2)) - {//Bragg - m_HXAi_bee[0] = 0; //ARF03062019 Change definition for the calculation of bee - m_HXAi_bee[1] = cos(thBrd) / m_dA; //ARF03062019 Change definition for the calculation of bee - m_HXAi_bee[2] = -sin(thBrd) / m_dA; //ARF03062019 Change definition for the calculation of bee - } - else if((m_uc == 3) || (m_uc == 4)) - {//Laue - m_HXAi_bee[0] = 0; //ARF03062019 Change definition for the calculation of bee - m_HXAi_bee[1] = -sin(thBrd) / m_dA; //IP14062019 for Laue pi/2 rotation of rec. vector from Bragg - m_HXAi_bee[2] = -cos(thBrd) / m_dA; //IP14062019 - } - // Conversion of polarization vector components to crystal frame: //e1X[0] = RLabXt[0][0]; //e1X[1] = RLabXt[1][0]; @@ -295,7 +279,7 @@ class srTOptCryst : public srTGenOptElem { **/ } - void FindDefOutFrameVect(double phEn, double hn, double ht, TVector3d& tv, TVector3d& sv, TVector3d& vZ1, TVector3d& vX1) + void FindDefOutFrameVect(double phEn, double hn, double ht, TVector3d& nv, TVector3d& tv, TVector3d& sv, TVector3d& vZ1, TVector3d& vX1) //ARF06062019 Introduce Vector3d& nv, { const double eAconv = 12398.4193009; double lam0Br = eAconv/phEn; //wavelength in [A] @@ -314,7 +298,7 @@ class srTOptCryst : public srTGenOptElem { if (absX1p > tolX1BackRefl) vX1p.Normalize(); else { vX1p.x = 1.; vX1p.y = 0.; vX1p.z = 0.; } - TVector3d nv = tv^sv; + //TVector3d nv = tv^sv; TVector3d vSt0_Rxl(sv.x, nv.x, tv.x); TVector3d vSt1_Rxl(sv.y, nv.y, tv.y); TVector3d vSt2_Rxl(sv.z, nv.z, tv.z); @@ -741,7 +725,7 @@ class srTOptCryst : public srTGenOptElem { TVector3d vZ1, vX1; //double hn = m_HXAi[1]; //OC180314: fix after conversation with A. Suvorov //double ht = m_HXAi[2]; - FindDefOutFrameVect(pRadAccessData->avgPhotEn, m_HXAi[1], m_HXAi[2], m_tv, m_sv, vZ1, vX1); + FindDefOutFrameVect(pRadAccessData->avgPhotEn, m_HXAi[1], m_HXAi[2], v_nv, m_tv, m_sv, vZ1, vX1); //Introduce v_nv FindRXtRef(vZ1, vX1, m_RLabXt, m_RXtRef); } @@ -754,8 +738,6 @@ class srTOptCryst : public srTGenOptElem { // return PropagateRadiationSingleE_Meth_0(pRadAccessData, 0); //} - //int PropagateRadiationSingleE_Meth_0(srTSRWRadStructAccessData* pRadAccessData, srTSRWRadStructAccessData* pPrevRadAccessData, void* pBuf=0) //OC06092019 - //OC01102019 (restored) int PropagateRadiationSingleE_Meth_0(srTSRWRadStructAccessData* pRadAccessData, srTSRWRadStructAccessData* pPrevRadAccessData) {//It works for many photon energies too (as in the case of Drift) //The "in-place" processing involving FFT for many photon energies greatly improves efficiency of the code for Time-/Frequency-Dependent simulations for FEL and pulsed lasers. @@ -828,8 +810,7 @@ class srTOptCryst : public srTGenOptElem { return 0; } - void RadPointModifier(srTEXZ& EXZ, srTEFieldPtrs& EPtrs, void* pBufVars=0) //OC29082019 - //void RadPointModifier(srTEXZ& EXZ, srTEFieldPtrs& EPtrs) + void RadPointModifier(srTEXZ& EXZ, srTEFieldPtrs& EPtrs) //void RadPointModifier_AngRepres(srTEXZ& EXZ, srTEFieldPtrs& EPtrs) {// E in eV; Length in m !!! // Operates on Angles side !!! @@ -927,14 +908,35 @@ class srTOptCryst : public srTGenOptElem { // Calculate direction cosine gamma0, reflection asymmetry parameter bee, // deviation parameter Adev and normalized deviation parameter zeeC (as // defined by Zachariasen gamma0, b, alpha and z.) - double gamma0 = -u0X[1]; + double gamma0; + + if(m_ug == 1) + { + gamma0 = -u0X[1]; + } + else if (m_ug == 2) + { + gamma0 = -u0X[1]; + } + double bee = 1. / (1. + m_HXAi[1] / k0wXAi[1]); - - //OC02102019 stopped updating here. - - //ARF03062019 suggestion: - //double bee = 1./(1. + (m_HXAi_bee[1]*m_RXtLab[1][1] + m_HXAi_bee[2]*m_RXtLab[2][1])/k0wXAi[1]); //ARF03062019 change from: double bee = 1. / (1. + m_HXAi[1] / k0wXAi[1]); + //ARF23052019 Change definition of bee + //double gamma0; + //double gammaH; + //double bee; + + //Bragg + //gamma0 = sin(pi/2 + thBrd + alphrd); + //gammaH = sin(-pi/2 + thBrd - alphrd); + //Laue + //if ((m_uc == 3) || (m_uc == 4)) + //{ + // gamma0 = sin(pi/2 - thBrd + alphrd); + // gammaH = sin(pi/2 + thBrd - alphrd); + //} + //bee = gamma0/gammaH; //ARF23052019 + double dotkH = k0wXAi[0] * m_HXAi[0] + k0wXAi[1] * m_HXAi[1] + k0wXAi[2] * m_HXAi[2]; double Adev = (2.*dotkH + 1. / (m_dA*m_dA)) / (k0Ai*k0Ai); complex zeeC = 0.5*((1. - bee)*m_psi0c + bee*Adev); @@ -948,60 +950,51 @@ class srTOptCryst : public srTGenOptElem { complex x2C = (-zeeC - sqrqzC) / m_psimhc; //complex ph1C = 2.*pi*iC*k0Ai*del1C*(m_thicum*1.E+04)/gamma0; //complex ph2C = 2.*pi*iC*k0Ai*del2C*(m_thicum*1.E+04)/gamma0; - complex aux_phC = 2.*pi*iC*k0Ai*(m_thicum*1.E+04) / gamma0; - complex ph1C = aux_phC*del1C; - complex ph2C = aux_phC*del2C; - - complex DHsgC, Cph1C, Cph2C, D0trsC; - if(real(ph1C) > logDBMax) DHsgC = x2C; - else if(real(ph2C) > logDBMax) DHsgC = x1C; - else + complex aux_phC = 2. * pi * iC * k0Ai * (m_thicum * 1.E+04) / gamma0; //m_thicum goes to A + complex ph1C = aux_phC * del1C; + complex ph2C = aux_phC * del2C; + complex Cph1C, Cph2C; + Cph1C = exp(ph1C); + Cph2C = exp(ph2C); + + complex DHsgC, D0trsC; + + if(m_ug == 1) //Bragg cases ARF080519 { - Cph1C = exp(ph1C); - Cph2C = exp(ph2C); - DHsgC = x1C*x2C*(Cph2C - Cph1C)/(Cph2C*x2C - Cph1C*x1C); - } - - //double re_ph1C = real(ph1C), re_ph2C = real(ph2C); - //if(re_ph1C > logDBMax) DHsgC = x2C; - //else if(re_ph2C > logDBMax) DHsgC = x1C; - //else - //{ - // if(re_ph1C < -logDBMax) - // { - // Cph1C = complex(0, 0); - // } - // else Cph1C = exp(ph1C); - - // if(re_ph2C < -logDBMax) - // { - // Cph2C = complex(0, 0); - // } - // else Cph2C = exp(ph2C); - - // if((re_ph1C < -logDBMax) && (re_ph2C < -logDBMax)) - // { - // DHsgC = complex(0, 0); //? - // } - // else DHsgC = x1C*x2C*(Cph2C - Cph1C)/(Cph2C*x2C - Cph1C*x1C); - //} - - if(m_itrans == 1) - { - // calculate the complex reflectivity D0trsC of the transmitted beam. - if(real(ph1C) > logDBMax) + //if(real(ph1C) > logDBMax) DHsgC = x2C; + //else if(real(ph2C) > logDBMax) DHsgC = x1C; + //else + //{ + // DHsgC = x1C * x2C * (Cph2C - Cph1C)/(Cph2C * x2C - Cph1C * x1C); + //} + DHsgC = x1C * x2C * (Cph2C - Cph1C)/(Cph2C * x2C - Cph1C * x1C); + + if(m_itrans == 1) { - Cph2C = exp(ph2C); - D0trsC = -Cph2C*(x2C - x1C)/x1C; + // calculate the complex reflectivity D0trsC of the transmitted beam. + //if(real(ph1C) > logDBMax) + //{ + // D0trsC = -Cph2C * (x2C - x1C)/x1C; + //} + //else if(real(ph2C) > logDBMax) + //{ + // D0trsC = +Cph1C * (x2C - x1C)/x2C; + //} + //else D0trsC = Cph1C * Cph2C * (x2C - x1C)/(Cph2C * x2C - Cph1C * x1C); + D0trsC = Cph1C * Cph2C * (x2C - x1C)/(Cph2C * x2C - Cph1C * x1C); } - else if(real(ph2C) > logDBMax) + } + //ARF 080519 See Zachariansen + else if(m_ug == 2) //Laue cases ARF080519 + { + DHsgC = x1C * x2C * (Cph1C - Cph2C)/(x2C - x1C); //Laue Reflection + // Transmission + if(m_itrans == 1) { - Cph1C = exp(ph1C); - D0trsC = +Cph1C*(x2C - x1C)/x2C; + D0trsC = (Cph1C * x2C - Cph2C * x1C)/(x2C - x1C); // Laue Transmission } - else D0trsC = Cph1C*Cph2C*(x2C - x1C)/(Cph2C*x2C - Cph1C*x1C); } - + // Calculate the complex reflectivity DHpiC for pi polarization: //cos2t = cos(2.*thBrd); //OC: moved to members m_cos2t queC = bee*m_psihc*m_psimhc*m_cos2t*m_cos2t; @@ -1012,33 +1005,45 @@ class srTOptCryst : public srTGenOptElem { x2C = (-zeeC - sqrqzC) / (m_psimhc*m_cos2t); ph1C = 2.*pi*iC*k0Ai*del1C*(m_thicum*1.E+04) / gamma0; ph2C = 2.*pi*iC*k0Ai*del2C*(m_thicum*1.E+04) / gamma0; + Cph1C = exp(ph1C); + Cph2C = exp(ph2C); complex DHpiC, D0trpC; - if(real(ph1C) > logDBMax) DHpiC = x2C; - else if(real(ph2C) > logDBMax) DHpiC = x1C; - else - { - Cph1C = exp(ph1C); - Cph2C = exp(ph2C); - DHpiC = x1C*x2C*(Cph2C - Cph1C)/(Cph2C*x2C - Cph1C*x1C); - } - if(m_itrans == 1) + if (m_ug == 1) //Bragg cases if introduced ARF080519 { - // calculate the complex reflectivity D0trpC of the transmitted beam. - if (real(ph1C) > logDBMax) + //if(real(ph1C) > logDBMax) DHpiC = x2C; + //else if(real(ph2C) > logDBMax) DHpiC = x1C; + //else + //{ + // DHpiC = x1C*x2C*(Cph2C - Cph1C)/(Cph2C*x2C - Cph1C*x1C); + //} + DHpiC = x1C*x2C*(Cph2C - Cph1C)/(Cph2C*x2C - Cph1C*x1C); + if(m_itrans == 1) { - Cph2C = exp(ph2C); - D0trpC = -Cph2C*(x2C - x1C)/x1C; + // calculate the complex reflectivity D0trpC of the transmitted beam. + //if (real(ph1C) > logDBMax) + //{ + // D0trpC = -Cph2C*(x2C - x1C)/x1C; + //} + //else if (real(ph2C) > logDBMax) + //{ + // D0trpC = +Cph1C*(x2C - x1C)/x2C; + //} + //else D0trpC = Cph1C*Cph2C*(x2C - x1C)/(Cph2C*x2C - Cph1C*x1C); + D0trpC = Cph1C*Cph2C*(x2C - x1C)/(Cph2C*x2C - Cph1C*x1C); } - else if (real(ph2C) > logDBMax) + } + else if(m_ug == 2) //Laue cases ARF080519 + { + DHpiC = x1C * x2C * (Cph1C - Cph2C)/(x2C - x1C); //Laue Reflection + //Transmission + if(m_itrans == 1) { - Cph1C = exp(ph1C); - D0trpC = +Cph1C*(x2C - x1C)/x2C; + D0trpC = (Cph1C * x2C - Cph2C*x1C)/(x2C - x1C); // Laue Transmission } - else D0trpC = Cph1C*Cph2C*(x2C - x1C)/(Cph2C*x2C - Cph1C*x1C); } - + // Calculate the diffracted amplitudes: complex EHSPCs = DHsgC*EInSPs; complex EHSPCp = DHpiC*EInSPp; @@ -1127,6 +1132,34 @@ class srTOptCryst : public srTGenOptElem { *(EPtrs.pEzRe) = (float)(Ety.real()); *(EPtrs.pEzIm) = (float)(Ety.imag()); } + else if(m_uc == 3) //Laue Reflection LD + { + complex Ehx = sgH[0]*EHSPCs + piH[0]*EHSPCp; + complex Ehy = sgH[1]*EHSPCs + piH[1]*EHSPCp; + complex Ehz = sgH[2]*EHSPCs + piH[2]*EHSPCp; //Longitudinal component not used (?) + + //Transverse components of the output electric field in the frame of the output beam: + *(EPtrs.pExRe) = (float)(asymFact*Ehx.real()); + *(EPtrs.pExIm) = (float)(asymFact*Ehx.imag()); + *(EPtrs.pEzRe) = (float)(asymFact*Ehy.real()); + *(EPtrs.pEzIm) = (float)(asymFact*Ehy.imag()); + } + else if(m_uc == 4) //Laue Transmission (i.e. Forward Bragg Diffraction (FLD)) + { + //DEBUG + //complex DHsgC_p_D0trsC = DHsgC + D0trsC; + //complex DHsgCae2_p_D0trsCae2 = DHsgC.real()*DHsgC.real() + DHsgC.imag()*DHsgC.imag() + D0trsC.real()*D0trsC.real() + D0trsC.imag()*D0trsC.imag(); + //END DEBUG + + complex Etx = m_InvPolTrn[0][0]*E0tSPs + m_InvPolTrn[0][1]*E0tSPp; + complex Ety = m_InvPolTrn[1][0]*E0tSPs + m_InvPolTrn[1][1]*E0tSPp; + + //Transverse components of the output electric field in the frame of the output beam: + *(EPtrs.pExRe) = (float)(Etx.real()); + *(EPtrs.pExIm) = (float)(Etx.imag()); + *(EPtrs.pEzRe) = (float)(Ety.real()); + *(EPtrs.pEzIm) = (float)(Ety.imag()); + } //OCTEST111014 //double testI1 = (*(EPtrs.pExRe))*(*(EPtrs.pExRe)) + (*(EPtrs.pExIm))*(*(EPtrs.pExIm)) +