LCOV - code coverage report
Current view: top level - EnergyPlus - SolarReflectionManager.cc (source / functions) Hit Total Coverage
Test: lcov.output.filtered Lines: 547 599 91.3 %
Date: 2023-01-17 19:17:23 Functions: 8 8 100.0 %

          Line data    Source code
       1             : // EnergyPlus, Copyright (c) 1996-2023, The Board of Trustees of the University of Illinois,
       2             : // The Regents of the University of California, through Lawrence Berkeley National Laboratory
       3             : // (subject to receipt of any required approvals from the U.S. Dept. of Energy), Oak Ridge
       4             : // National Laboratory, managed by UT-Battelle, Alliance for Sustainable Energy, LLC, and other
       5             : // contributors. All rights reserved.
       6             : //
       7             : // NOTICE: This Software was developed under funding from the U.S. Department of Energy and the
       8             : // U.S. Government consequently retains certain rights. As such, the U.S. Government has been
       9             : // granted for itself and others acting on its behalf a paid-up, nonexclusive, irrevocable,
      10             : // worldwide license in the Software to reproduce, distribute copies to the public, prepare
      11             : // derivative works, and perform publicly and display publicly, and to permit others to do so.
      12             : //
      13             : // Redistribution and use in source and binary forms, with or without modification, are permitted
      14             : // provided that the following conditions are met:
      15             : //
      16             : // (1) Redistributions of source code must retain the above copyright notice, this list of
      17             : //     conditions and the following disclaimer.
      18             : //
      19             : // (2) Redistributions in binary form must reproduce the above copyright notice, this list of
      20             : //     conditions and the following disclaimer in the documentation and/or other materials
      21             : //     provided with the distribution.
      22             : //
      23             : // (3) Neither the name of the University of California, Lawrence Berkeley National Laboratory,
      24             : //     the University of Illinois, U.S. Dept. of Energy nor the names of its contributors may be
      25             : //     used to endorse or promote products derived from this software without specific prior
      26             : //     written permission.
      27             : //
      28             : // (4) Use of EnergyPlus(TM) Name. If Licensee (i) distributes the software in stand-alone form
      29             : //     without changes from the version obtained under this License, or (ii) Licensee makes a
      30             : //     reference solely to the software portion of its product, Licensee must refer to the
      31             : //     software as "EnergyPlus version X" software, where "X" is the version number Licensee
      32             : //     obtained under this License and may not use a different name for the software. Except as
      33             : //     specifically required in this Section (4), Licensee shall not use in a company name, a
      34             : //     product name, in advertising, publicity, or other promotional activities any name, trade
      35             : //     name, trademark, logo, or other designation of "EnergyPlus", "E+", "e+" or confusingly
      36             : //     similar designation, without the U.S. Department of Energy's prior written consent.
      37             : //
      38             : // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
      39             : // IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY
      40             : // AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
      41             : // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
      42             : // CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
      43             : // SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
      44             : // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
      45             : // OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
      46             : // POSSIBILITY OF SUCH DAMAGE.
      47             : 
      48             : // C++ Headers
      49             : #include <cassert>
      50             : #include <cmath>
      51             : 
      52             : // ObjexxFCL Headers
      53             : #include <ObjexxFCL/Fmath.hh>
      54             : 
      55             : // EnergyPlus Headers
      56             : #include <EnergyPlus/Construction.hh>
      57             : #include <EnergyPlus/Data/EnergyPlusData.hh>
      58             : #include <EnergyPlus/DataEnvironment.hh>
      59             : #include <EnergyPlus/DataHeatBalance.hh>
      60             : #include <EnergyPlus/DataSurfaces.hh>
      61             : #include <EnergyPlus/DataSystemVariables.hh>
      62             : #include <EnergyPlus/DataVectorTypes.hh>
      63             : #include <EnergyPlus/DisplayRoutines.hh>
      64             : #include <EnergyPlus/General.hh>
      65             : #include <EnergyPlus/PierceSurface.hh>
      66             : #include <EnergyPlus/ScheduleManager.hh>
      67             : #include <EnergyPlus/SolarReflectionManager.hh>
      68             : #include <EnergyPlus/SolarShading.hh>
      69             : 
      70             : namespace EnergyPlus {
      71             : 
      72             : namespace SolarReflectionManager {
      73             : 
      74             :     // MODULE INFORMATION
      75             :     //       AUTHOR         Fred Winkelmann
      76             :     //       DATE WRITTEN   September 2003
      77             :     //       MODIFIED       May 2004, FCW: modify calculation of receiving point location on a
      78             :     //                        receiving surface so can handle surface of any number of vertices
      79             :     //                        (previously restricted to 3- or 4-sided surfaces).
      80             :     //       RE-ENGINEERED  na
      81             : 
      82             :     // PURPOSE OF THIS MODULE:
      83             :     // Manages the calculation of factors for solar reflected from obstructions and ground.
      84             : 
      85             :     // METHODOLOGY EMPLOYED:
      86             :     // REFERENCES: na
      87             : 
      88             :     // OTHER NOTES: na
      89             : 
      90             :     // Using/Aliasing
      91             :     using namespace DataHeatBalance;
      92             :     using namespace DataSurfaces;
      93             :     using namespace ScheduleManager;
      94             :     using namespace DataEnvironment;
      95             : 
      96             :     using namespace DataVectorTypes;
      97             : 
      98           9 :     void InitSolReflRecSurf(EnergyPlusData &state)
      99             :     {
     100             : 
     101             :         // SUBROUTINE INFORMATION:
     102             :         //       AUTHOR         Fred Winkelmann
     103             :         //       DATE WRITTEN   September 2003
     104             :         //       MODIFIED       na
     105             :         //       RE-ENGINEERED  na
     106             : 
     107             :         // PURPOSE OF THIS SUBROUTINE:
     108             :         // Initializes the derived type SolReflRecSurf, which contains information
     109             :         // needed to calculate factors for solar reflection from obstructions and ground.
     110             : 
     111             :         // SUBROUTINE LOCAL VARIABLE DECLARATIONS:
     112             :         int SurfNum;            // Surface number
     113             :         int RecSurfNum;         // Receiving surface number
     114             :         int loop;               // DO loop indices
     115             :         int loop1;              // DO loop indices
     116             :         int loopA;              // DO loop indices
     117             :         int loopB;              // DO loop indices
     118             :         int ObsSurfNum;         // Surface number of an obstruction
     119             :         bool ObsBehindRec;      // True if an obstruction is entirely behind a receiving surface
     120             :         bool ObsHasView;        // True if view between receiving surface and heat trans surf obstruction
     121          18 :         Vector3<Real64> RecVec; // First vertex of a receiving surface (m)
     122          18 :         Vector3<Real64> ObsVec; // A vertex of a candidate obstructing surface (m)
     123          18 :         Vector3<Real64> VecAB;  // Vector from receiving surface vertex to obstruction surface vertex (m)
     124          18 :         Vector3<Real64> HitPt;  // Hit point (m)
     125             :         Real64 DotProd;         // Dot product of vectors (m2)
     126             :         int RecPtNum;           // Receiving point number
     127             :         // unused  REAL(r64)         :: SumX                 ! Sum of X (or Y or Z) coordinate values of a surface
     128             :         // unused  REAL(r64)         :: SumY                 ! Sum of X (or Y or Z) coordinate values of a surface
     129             :         // unused  REAL(r64)         :: SumZ                 ! Sum of X (or Y or Z) coordinate values of a surface
     130             :         Real64 PhiSurf;   // Altitude of normal to receiving surface (radians)
     131             :         Real64 ThetaSurf; // Azimuth of normal to receiving surface (radians)
     132             :         Real64 PhiMin;    // Minimum and maximum values of ray altitude angle (radians)
     133             :         Real64 PhiMax;    // Minimum and maximum values of ray altitude angle (radians)
     134             :         Real64 ThetaMin;  // Minimum and maximum values of ray azimuth angle (radians)
     135             :         Real64 ThetaMax;  // Minimum and maximum values of ray azimuth angle (radians)
     136             :         Real64 Phi;       // Ray altitude angle, increment, sine, and cosine
     137             :         Real64 DPhi;      // Ray altitude angle, increment, sine, and cosine
     138             :         Real64 SPhi;      // Ray altitude angle, increment, sine, and cosine
     139             :         Real64 CPhi;      // Ray altitude angle, increment, sine, and cosine
     140             :         Real64 Theta;     // Ray azimuth angle and increment
     141             :         Real64 DTheta;    // Ray azimuth angle and increment
     142             :         int IPhi;         // Ray altitude angle and azimuth angle indices
     143             :         int ITheta;       // Ray altitude angle and azimuth angle indices
     144             :         // unused  REAL(r64)         :: APhi                 ! Intermediate variable
     145             :         int RayNum;                   // Ray number
     146          18 :         Vector3<Real64> URay;         // Unit vector along ray pointing away from receiving surface
     147             :         Real64 CosIncAngRay;          // Cosine of angle of incidence of ray on receiving surface
     148             :         Real64 dOmega;                // Solid angle associated with a ray
     149             :         bool hit;                     // True iff obstruction is hit
     150             :         int TotObstructionsHit;       // Number of obstructions hit by a ray
     151             :         Real64 HitDistance;           // Distance from receiving point to hit point for a ray (m)
     152             :         int NearestHitSurfNum;        // Surface number of nearest obstruction hit by a ray
     153          18 :         Vector3<Real64> NearestHitPt; // Nearest hit pit for a ray (m)
     154             :         Real64 NearestHitDistance;    // Distance from receiving point to nearest hit point for a ray (m)
     155             :         int ObsSurfNumToSkip;         // Surface number of obstruction to be ignored
     156          18 :         Vector3<Real64> RecPt;        // Receiving point (m)
     157          18 :         Vector3<Real64> RayVec;       // Unit vector along ray
     158          18 :         Vector3<Real64> Vec1;         // Vectors between hit surface vertices (m)
     159          18 :         Vector3<Real64> Vec2;         // Vectors between hit surface vertices (m)
     160          18 :         Vector3<Real64> VNorm;        // For a hit surface, unit normal vector pointing into the hemisphere
     161             :         // containing the receiving point
     162             :         int ObsConstrNum; // Construction number of obstruction; = 0 if a shading surface
     163             :         Real64 Alfa;      // Direction angles for ray heading towards the ground (radians)
     164             :         Real64 Beta;
     165             :         Real64 HorDis;               // Distance between ground hit point and proj'n of receiving pt onto ground (m)
     166          18 :         Vector3<Real64> GroundHitPt; // Coordinates of ground hit point
     167             :         // unused  REAL(r64)         :: ArgASin
     168             :         Real64 ACosTanTan;
     169             :         int J;           // DO loop indices
     170             :         int K;           // DO loop indices
     171             :         int NumRecPts;   // Number of surface receiving points for reflected solar radiation
     172             :         Real64 VertexWt; // Vertex weighting factor for calculating receiving points
     173             : 
     174           9 :         static Vector3<Real64> const unit_z(0.0, 0.0, 1.0);
     175           9 :         static Vector3<Real64> const zero3(0.0);
     176             : 
     177             :         // Find number of surfaces that are sun-exposed exterior building heat transfer surfaces.
     178             :         // These are candidates for receiving solar reflected from obstructions and ground.
     179             :         // CR 7640.  12/3/2008 BG simplified logic to allow for Other Side Conditions Modeled boundary condition.
     180             :         //           and solar collectors on shading surfaces that need this.
     181             : 
     182             :         // shading surfaces have ExtSolar = False, so they are not included in TotSolReflRecSurf
     183           9 :         state.dataSolarReflectionManager->TotSolReflRecSurf = 0;
     184         101 :         for (SurfNum = 1; SurfNum <= state.dataSurface->TotSurfaces; ++SurfNum) {
     185          92 :             if (state.dataSurface->Surface(SurfNum).ExtSolar) {
     186          47 :                 ++state.dataSolarReflectionManager->TotSolReflRecSurf;
     187             :             }
     188             :         }
     189             : 
     190             :         // TH 3/29/2010. ShadowSurfPossibleReflector is not used!
     191             :         // Set flag that determines whether a surface can be an exterior reflector
     192             :         // DO SurfNum = 1,TotSurfaces
     193             :         //  Surface(SurfNum)%ShadowSurfPossibleReflector = .FALSE.
     194             :         // Exclude non-exterior heat transfer surfaces (but not OtherSideCondModeledExt = -4 CR7640)
     195             :         //  IF(Surface(SurfNum)%HeatTransSurf .AND. Surface(SurfNum)%ExtBoundCond > 0 ) CYCLE
     196             :         //  IF(Surface(SurfNum)%HeatTransSurf .AND. Surface(SurfNum)%ExtBoundCond == Ground) CYCLE
     197             :         //  IF(Surface(SurfNum)%HeatTransSurf .AND. Surface(SurfNum)%ExtBoundCond == OtherSideCoefNoCalcExt) CYCLE
     198             :         //  IF(Surface(SurfNum)%HeatTransSurf .AND. Surface(SurfNum)%ExtBoundCond == OtherSideCoefCalcExt) CYCLE
     199             : 
     200             :         // Exclude daylighting shelves. A separate solar reflection calculation is done for these.
     201             :         //  IF(Surface(SurfNum)%Shelf > 0) CYCLE
     202             : 
     203             :         // Exclude duplicate shading surfaces
     204             :         // TH 3/24/2010. Why? a mirror shading surface can reflect solar (either beam or diffuse)
     205             :         //  can use a flag like Surface(SurfNum)%Mirrored (True or False) to avoid string comparison
     206             :         //   and to allow surface names starting with 'Mir'
     207             :         // IF(Surface(SurfNum)%Name(1:3) == 'Mir') CYCLE
     208             :         //  IF(Surface(SurfNum)%MirroredSurf) CYCLE
     209             : 
     210             :         //  Surface(SurfNum)%ShadowSurfPossibleReflector = .TRUE.
     211             :         // END DO
     212             : 
     213           9 :         if (state.dataSolarReflectionManager->TotSolReflRecSurf == 0) {
     214           0 :             ShowWarningError(state, "Calculation of solar reflected from obstructions has been requested but there");
     215           0 :             ShowContinueError(state, "are no building surfaces that can receive reflected solar. Calculation will not be done.");
     216           0 :             state.dataSurface->CalcSolRefl = false;
     217           0 :             return;
     218             :         }
     219             : 
     220             :         // Should this be moved up front?
     221           9 :         if (state.dataEnvrn->IgnoreSolarRadiation) {
     222           0 :             state.dataSolarReflectionManager->TotSolReflRecSurf = 0;
     223           0 :             state.dataSurface->CalcSolRefl = false;
     224           0 :             return;
     225             :         }
     226             : 
     227           9 :         state.dataSolarReflectionManager->SolReflRecSurf.allocate(state.dataSolarReflectionManager->TotSolReflRecSurf);
     228             : 
     229           9 :         state.dataSurface->SurfReflFacBmToDiffSolObs.dimension(24, state.dataSurface->TotSurfaces, 0.0);
     230           9 :         state.dataSurface->SurfReflFacBmToDiffSolGnd.dimension(24, state.dataSurface->TotSurfaces, 0.0);
     231           9 :         state.dataSurface->SurfReflFacBmToBmSolObs.dimension(24, state.dataSurface->TotSurfaces, 0.0);
     232           9 :         state.dataSurface->SurfReflFacSkySolObs.dimension(state.dataSurface->TotSurfaces, 0.0);
     233           9 :         state.dataSurface->SurfReflFacSkySolGnd.dimension(state.dataSurface->TotSurfaces, 0.0);
     234           9 :         state.dataSurface->SurfCosIncAveBmToBmSolObs.dimension(24, state.dataSurface->TotSurfaces, 0.0);
     235             : 
     236             :         // Only surfaces with sun exposure can receive solar reflection from ground or onstructions
     237             :         //  Shading surfaces are always not exposed to solar (ExtSolar = False)
     238           9 :         RecSurfNum = 0;
     239         101 :         for (SurfNum = 1; SurfNum <= state.dataSurface->TotSurfaces; ++SurfNum) {
     240          92 :             state.dataSurface->SurfShadowRecSurfNum(SurfNum) = 0;
     241          92 :             if (state.dataSurface->Surface(SurfNum).ExtSolar) {
     242          47 :                 ++RecSurfNum;
     243          47 :                 state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).SurfNum = SurfNum;
     244          47 :                 state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).SurfName = state.dataSurface->Surface(SurfNum).Name;
     245          47 :                 state.dataSurface->SurfShadowRecSurfNum(SurfNum) = RecSurfNum;
     246             : 
     247             :                 // Warning if any receiving surface vertex is below ground level, taken to be at Z = 0 in absolute coords
     248         235 :                 for (loop = 1; loop <= state.dataSurface->Surface(SurfNum).Sides; ++loop) {
     249         188 :                     if (state.dataSurface->Surface(SurfNum).Vertex(loop).z < state.dataSurface->GroundLevelZ) {
     250           0 :                         ShowWarningError(
     251           0 :                             state, "Calculation of reflected solar onto surface=" + state.dataSurface->Surface(SurfNum).Name + " may be inaccurate");
     252           0 :                         ShowContinueError(state, "because it has one or more vertices below ground level.");
     253           0 :                         break;
     254             :                     }
     255             :                 }
     256             :             }
     257             :         }
     258             : 
     259             :         // Get MaxRecPts for allocating SolReflRecSurf arrays that depend on number of receiving points
     260           9 :         state.dataSurface->MaxRecPts = 1;
     261          56 :         for (RecSurfNum = 1; RecSurfNum <= state.dataSolarReflectionManager->TotSolReflRecSurf; ++RecSurfNum) {
     262          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumRecPts =
     263          47 :                 state.dataSurface->Surface(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).SurfNum).Sides;
     264          47 :             if (state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumRecPts > state.dataSurface->MaxRecPts)
     265           9 :                 state.dataSurface->MaxRecPts = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumRecPts;
     266             :         }
     267             : 
     268           9 :         state.dataSurface->MaxReflRays = AltAngStepsForSolReflCalc * AzimAngStepsForSolReflCalc;
     269          56 :         for (RecSurfNum = 1; RecSurfNum <= state.dataSolarReflectionManager->TotSolReflRecSurf; ++RecSurfNum) {
     270          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec = 0.0;
     271          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPt.dimension(state.dataSurface->MaxRecPts, zero3);
     272          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RayVec.dimension(state.dataSurface->MaxReflRays, zero3);
     273          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).CosIncAngRay.dimension(state.dataSurface->MaxReflRays, 0.0);
     274          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).dOmegaRay.dimension(state.dataSurface->MaxReflRays, 0.0);
     275          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum)
     276          47 :                 .HitPt.dimension(state.dataSurface->MaxReflRays, state.dataSurface->MaxRecPts, zero3);
     277          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum)
     278          47 :                 .HitPtSurfNum.dimension(state.dataSurface->MaxReflRays, state.dataSurface->MaxRecPts, 0.0);
     279          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum)
     280          47 :                 .HitPtSolRefl.dimension(state.dataSurface->MaxReflRays, state.dataSurface->MaxRecPts, 0.0);
     281          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum)
     282          47 :                 .RecPtHitPtDis.dimension(state.dataSurface->MaxReflRays, state.dataSurface->MaxRecPts, 0.0);
     283          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum)
     284          47 :                 .HitPtNormVec.dimension(state.dataSurface->MaxReflRays, state.dataSurface->MaxRecPts, zero3);
     285          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).PossibleObsSurfNums.dimension(state.dataSurface->TotSurfaces, 0);
     286             :         }
     287             : 
     288          56 :         for (RecSurfNum = 1; RecSurfNum <= state.dataSolarReflectionManager->TotSolReflRecSurf; ++RecSurfNum) {
     289          47 :             SurfNum = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).SurfNum;
     290             :             // Outward norm to receiving surface
     291          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec = state.dataSurface->Surface(SurfNum).OutNormVec;
     292          47 :             RecVec = state.dataSurface->Surface(SurfNum).Vertex(1);
     293             :             // Loop over all surfaces and find those that can be obstructing surfaces for this receiving surf
     294          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumPossibleObs = 0;
     295         683 :             for (ObsSurfNum = 1; ObsSurfNum <= state.dataSurface->TotSurfaces; ++ObsSurfNum) {
     296             :                 // Exclude the receiving surface itself and its base surface (if it has one)
     297         636 :                 if (ObsSurfNum == SurfNum || ObsSurfNum == state.dataSurface->Surface(SurfNum).BaseSurf) continue;
     298             :                 // Exclude non-exterior heat transfer surfaces
     299         574 :                 if (state.dataSurface->Surface(ObsSurfNum).HeatTransSurf && state.dataSurface->Surface(ObsSurfNum).ExtBoundCond != 0) continue;
     300             :                 // Exclude duplicate shading surfaces
     301             :                 // IF(Surface(ObsSurfNum)%Name(1:3) == 'Mir') CYCLE
     302             :                 // TH2 CR8959
     303             :                 // IF(Surface(ObsSurfNum)%MirroredSurf) CYCLE
     304             : 
     305             :                 // Exclude surfaces that are entirely behind the receiving surface.This is true if dot products of the
     306             :                 // rec. surface outward normal and vector from first vertex of rec. surface and each vertex of
     307             :                 // obstructing surface are all negative.
     308         353 :                 ObsBehindRec = true;
     309        1611 :                 for (loop = 1; loop <= state.dataSurface->Surface(ObsSurfNum).Sides; ++loop) {
     310        1303 :                     ObsVec = state.dataSurface->Surface(ObsSurfNum).Vertex(loop);
     311        1303 :                     DotProd = dot(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec, ObsVec - RecVec);
     312             :                     // CR8251      IF(DotProd > 0.01d0) THEN  ! This obstructing-surface vertex is not behind receiving surface
     313        1303 :                     if (DotProd > 1.0e-6) { // This obstructing-surface vertex is not behind receiving surface
     314          45 :                         ObsBehindRec = false;
     315          45 :                         break;
     316             :                     }
     317             :                 }
     318         353 :                 if (ObsBehindRec) continue;
     319             : 
     320             :                 // Exclude heat transfer surfaces that have no view with the receiving surface.
     321             :                 // There is view if: for at least one vector VecAB from a receiving surface vertex to
     322             :                 // a vertex of a potential obstructing surface that satisfies VecAB.nA > 0.0 and VecAB.nB < 0.0,
     323             :                 // where nA and nB are the outward normal to the receiving and obstructing surface, resp.
     324          45 :                 if (state.dataSurface->Surface(ObsSurfNum).HeatTransSurf) {
     325          17 :                     ObsHasView = false;
     326          67 :                     for (loopA = 1; loopA <= state.dataSurface->Surface(SurfNum).Sides; ++loopA) {
     327         260 :                         for (loopB = 1; loopB <= state.dataSurface->Surface(ObsSurfNum).Sides; ++loopB) {
     328         210 :                             VecAB = (state.dataSurface->Surface(ObsSurfNum).Vertex(loopB) - state.dataSurface->Surface(SurfNum).Vertex(loopA));
     329         360 :                             if (dot(VecAB, state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec) > 0.0 &&
     330         150 :                                 dot(VecAB, state.dataSurface->Surface(ObsSurfNum).OutNormVec) < 0.0) {
     331           6 :                                 ObsHasView = true;
     332           6 :                                 break;
     333             :                             }
     334             :                         }
     335          56 :                         if (ObsHasView) break;
     336             :                     }
     337          17 :                     if (!ObsHasView) continue;
     338             :                 }
     339             : 
     340             :                 // This is a possible obstructing surface for this receiving surface
     341          34 :                 ++state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumPossibleObs;
     342          34 :                 state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum)
     343          68 :                     .PossibleObsSurfNums(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumPossibleObs) = ObsSurfNum;
     344             :             }
     345             : 
     346             :             // Get coordinates of receiving points on this receiving surface. The number of receiving points
     347             :             // is equal to the number of surface vertices (3 or higher).
     348             : 
     349          47 :             NumRecPts = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumRecPts;
     350         235 :             for (J = 1; J <= NumRecPts; ++J) {
     351         188 :                 state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPt(J) = 0.0;
     352         940 :                 for (K = 1; K <= NumRecPts; ++K) {
     353         752 :                     if (NumRecPts == 3) { // Receiving surface is a triangle
     354           0 :                         VertexWt = 0.2;
     355           0 :                         if (K == J) VertexWt = 0.6;
     356             :                     } else { // Receiving surface has 4 or more vertices
     357         752 :                         VertexWt = 1.0 / (2.0 * NumRecPts);
     358         752 :                         if (K == J) VertexWt = (NumRecPts + 1.0) / (2.0 * NumRecPts);
     359             :                     }
     360         752 :                     state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPt(J).x +=
     361         752 :                         VertexWt * state.dataSurface->Surface(SurfNum).Vertex(K).x;
     362         752 :                     state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPt(J).y +=
     363         752 :                         VertexWt * state.dataSurface->Surface(SurfNum).Vertex(K).y;
     364         752 :                     state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPt(J).z +=
     365         752 :                         VertexWt * state.dataSurface->Surface(SurfNum).Vertex(K).z;
     366             :                 }
     367             :             }
     368             : 
     369             :             // Create rays going outward from receiving surface. The same rays will be used at each receiving point.
     370             :             // The rays are used in calculating diffusely reflected solar incident on receiving surface.
     371             : 
     372             :             // Divide hemisphere around receiving surface into elements of altitude Phi and
     373             :             // azimuth Theta and create ray unit vector at each Phi,Theta pair in front of the surface.
     374             :             // Phi = 0 at the horizon; Phi = Pi/2 at the zenith
     375             : 
     376          47 :             PhiSurf = std::asin(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec.z);
     377          47 :             Real64 const tan_PhiSurf = std::tan(PhiSurf);
     378          47 :             Real64 const sin_PhiSurf = std::sin(PhiSurf);
     379          47 :             Real64 const cos_PhiSurf = std::cos(PhiSurf);
     380             : 
     381          72 :             if (std::abs(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec.x) > 1.0e-5 ||
     382          25 :                 std::abs(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec.y) > 1.0e-5) {
     383          38 :                 ThetaSurf = std::atan2(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec.y,
     384          38 :                                        state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec.x);
     385             :             } else {
     386           9 :                 ThetaSurf = 0.0;
     387             :             }
     388          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).PhiNormVec = PhiSurf;
     389          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).ThetaNormVec = ThetaSurf;
     390          47 :             PhiMin = max(-DataGlobalConstants::PiOvr2, PhiSurf - DataGlobalConstants::PiOvr2);
     391          47 :             PhiMax = min(DataGlobalConstants::PiOvr2, PhiSurf + DataGlobalConstants::PiOvr2);
     392          47 :             DPhi = (PhiMax - PhiMin) / AltAngStepsForSolReflCalc;
     393          47 :             RayNum = 0;
     394             : 
     395             :             // Altitude loop
     396         517 :             for (IPhi = 1; IPhi <= AltAngStepsForSolReflCalc; ++IPhi) {
     397         470 :                 Phi = PhiMin + (IPhi - 0.5) * DPhi;
     398         470 :                 SPhi = std::sin(Phi);
     399         470 :                 CPhi = std::cos(Phi);
     400             :                 // Third component of ray unit vector in (Theta,Phi) direction
     401         470 :                 URay(3) = SPhi;
     402             : 
     403         470 :                 if (PhiSurf >= 0.0) {
     404         450 :                     if (Phi >= DataGlobalConstants::PiOvr2 - PhiSurf) {
     405          70 :                         ThetaMin = -DataGlobalConstants::Pi;
     406          70 :                         ThetaMax = DataGlobalConstants::Pi;
     407             :                     } else {
     408         380 :                         ACosTanTan = std::acos(-std::tan(Phi) * tan_PhiSurf);
     409         380 :                         ThetaMin = ThetaSurf - std::abs(ACosTanTan);
     410         380 :                         ThetaMax = ThetaSurf + std::abs(ACosTanTan);
     411             :                     }
     412             : 
     413             :                 } else { // PhiSurf < 0.0
     414          20 :                     if (Phi <= -PhiSurf - DataGlobalConstants::PiOvr2) {
     415          20 :                         ThetaMin = -DataGlobalConstants::Pi;
     416          20 :                         ThetaMax = DataGlobalConstants::Pi;
     417             :                     } else {
     418           0 :                         ACosTanTan = std::acos(-std::tan(Phi) * tan_PhiSurf);
     419           0 :                         ThetaMin = ThetaSurf - std::abs(ACosTanTan);
     420           0 :                         ThetaMax = ThetaSurf + std::abs(ACosTanTan);
     421             :                     }
     422             :                 }
     423             : 
     424         470 :                 DTheta = (ThetaMax - ThetaMin) / AzimAngStepsForSolReflCalc;
     425         470 :                 dOmega = CPhi * DTheta * DPhi;
     426             : 
     427             :                 // Azimuth loop
     428        4700 :                 for (ITheta = 1; ITheta <= AzimAngStepsForSolReflCalc; ++ITheta) {
     429        4230 :                     Theta = ThetaMin + (ITheta - 0.5) * DTheta;
     430        4230 :                     URay.x = CPhi * std::cos(Theta);
     431        4230 :                     URay.y = CPhi * std::sin(Theta);
     432             :                     // Cosine of angle of incidence of ray on receiving surface
     433        4230 :                     CosIncAngRay = SPhi * sin_PhiSurf + CPhi * cos_PhiSurf * std::cos(Theta - ThetaSurf);
     434        4230 :                     if (CosIncAngRay < 0.0) continue; // Ray is behind receiving surface (although there shouldn't be any)
     435        4230 :                     ++RayNum;
     436        4230 :                     state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RayVec(RayNum) = URay;
     437        4230 :                     state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).CosIncAngRay(RayNum) = CosIncAngRay;
     438        4230 :                     state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).dOmegaRay(RayNum) = dOmega;
     439             :                 } // End of azimuth loop
     440             : 
     441             :             } // End of altitude loop
     442          47 :             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumReflRays = RayNum;
     443             : 
     444             :         } // End of loop over receiving surfaces
     445             : 
     446             :         // Loop again over receiving surfaces and, for each ray, get hit point and info associated with that point
     447             :         // (hit point = point that ray intersects nearest obstruction, or, if ray is downgoing and hits no
     448             :         // obstructions, point that ray intersects ground plane).
     449             : 
     450          56 :         for (RecSurfNum = 1; RecSurfNum <= state.dataSolarReflectionManager->TotSolReflRecSurf; ++RecSurfNum) {
     451          47 :             SurfNum = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).SurfNum;
     452         235 :             for (RecPtNum = 1; RecPtNum <= state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumRecPts; ++RecPtNum) {
     453         188 :                 RecPt = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPt(RecPtNum);
     454       17108 :                 for (RayNum = 1; RayNum <= state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumReflRays; ++RayNum) {
     455             :                     // Loop over possible obstructions. If ray hits one or more obstructions get hit point on closest obstruction.
     456             :                     // If ray hits no obstructions and is going upward set HitPointSurfNum = 0.
     457             :                     // If ray hits no obstructions and is going downward set HitPointSurfNum = -1 and get hit point on ground.
     458       16920 :                     TotObstructionsHit = 0;
     459       16920 :                     NearestHitSurfNum = 0;
     460       16920 :                     NearestHitDistance = 1.0e+8;
     461       16920 :                     ObsSurfNumToSkip = 0;
     462       16920 :                     RayVec = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RayVec(RayNum);
     463       29160 :                     for (loop1 = 1; loop1 <= state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumPossibleObs; ++loop1) {
     464             :                         // Surface number of this obstruction
     465       12240 :                         ObsSurfNum = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).PossibleObsSurfNums(loop1);
     466             :                         // If a window was hit previously (see below), ObsSurfNumToSkip was set to the window's base surface in order
     467             :                         // to remove that surface from consideration as a hit surface for this ray
     468       12240 :                         if (ObsSurfNum == ObsSurfNumToSkip) continue;
     469             :                         // Determine if this ray hits ObsSurfNum (in which case hit is true) and, if so, what the
     470             :                         // distance from the receiving point to the hit point is
     471             :                         PierceSurface(state, ObsSurfNum, RecPt, RayVec, HitPt, hit);
     472       12240 :                         if (hit) {
     473             :                             // added TH 3/29/2010 to set ObsSurfNumToSkip
     474        1432 :                             if (state.dataSurface->Surface(ObsSurfNum).Class == SurfaceClass::Window) {
     475         104 :                                 ObsSurfNumToSkip = state.dataSurface->Surface(ObsSurfNum).BaseSurf;
     476             :                             }
     477             : 
     478             :                             // If obstruction is a window and its base surface is the nearest obstruction hit so far,
     479             :                             // set NearestHitSurfNum to this window. Note that in this case NearestHitDistance has already
     480             :                             // been calculated, so does not have to be recalculated.
     481        1536 :                             if (state.dataSurface->Surface(ObsSurfNum).Class == SurfaceClass::Window &&
     482         104 :                                 state.dataSurface->Surface(ObsSurfNum).BaseSurf == NearestHitSurfNum) {
     483         104 :                                 NearestHitSurfNum = ObsSurfNum;
     484             :                             } else {
     485        1328 :                                 ++TotObstructionsHit;
     486             :                                 // Distance from receiving point to hit point
     487        1328 :                                 HitDistance = distance(HitPt, RecPt);
     488             :                                 // Reset NearestHitSurfNum and NearestHitDistance if this hit point is closer than previous closest
     489        1328 :                                 if (HitDistance < NearestHitDistance) {
     490         834 :                                     NearestHitDistance = HitDistance;
     491         834 :                                     NearestHitSurfNum = ObsSurfNum;
     492         834 :                                     NearestHitPt = HitPt;
     493         494 :                                 } else if (HitDistance == NearestHitDistance) { // TH2 CR8959
     494             :                                     // Ray hits mirrored surfaces. Choose the surface facing the ray.
     495         494 :                                     if (dot(state.dataSurface->Surface(ObsSurfNum).OutNormVec, RayVec) <= 0.0) {
     496         494 :                                         NearestHitSurfNum = ObsSurfNum;
     497             :                                     }
     498             :                                 }
     499             :                             }
     500             :                         } // End of check if obstruction was hit
     501             :                     }     // End of loop over possible obstructions for this ray
     502             : 
     503       16920 :                     if (TotObstructionsHit > 0) {
     504             :                         // One or more obstructions were hit by this ray
     505         824 :                         state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPtSurfNum(RayNum, RecPtNum) = NearestHitSurfNum;
     506         824 :                         state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPtHitPtDis(RayNum, RecPtNum) = NearestHitDistance;
     507         824 :                         state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPt(RayNum, RecPtNum) = NearestHitPt;
     508             :                         // For hit surface, calculate unit normal vector pointing into the hemisphere
     509             :                         // containing the receiving point
     510         824 :                         Vec1 = (state.dataSurface->Surface(NearestHitSurfNum).Vertex(1) - state.dataSurface->Surface(NearestHitSurfNum).Vertex(3));
     511         824 :                         Vec2 = (state.dataSurface->Surface(NearestHitSurfNum).Vertex(2) - state.dataSurface->Surface(NearestHitSurfNum).Vertex(3));
     512         824 :                         VNorm = cross(Vec1, Vec2);
     513         824 :                         VNorm.normalize(); // Do Handle magnitude==0
     514         824 :                         if (dot(VNorm, -RayVec) < 0.0) VNorm = -VNorm;
     515         824 :                         state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPtNormVec(RayNum, RecPtNum) = VNorm;
     516             :                         // Get solar and visible beam-to-diffuse reflectance at nearest hit point
     517         824 :                         ObsConstrNum = state.dataSurface->Surface(NearestHitSurfNum).Construction;
     518         824 :                         if (ObsConstrNum > 0) {
     519             :                             // Exterior building surface is nearest hit
     520         340 :                             if (!state.dataConstruction->Construct(ObsConstrNum).TypeIsWindow) {
     521             :                                 // Obstruction is not a window, i.e., is an opaque surface
     522         236 :                                 state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPtSolRefl(RayNum, RecPtNum) =
     523         236 :                                     1.0 - state.dataConstruction->Construct(ObsConstrNum).OutsideAbsorpSolar;
     524             :                             } else {
     525             :                                 // Obstruction is a window. Assume it is bare so that there is no beam-to-diffuse reflection
     526             :                                 // (beam-to-beam reflection is calculated in subroutine CalcBeamSolSpecularReflFactors).
     527         104 :                                 state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPtSolRefl(RayNum, RecPtNum) = 0.0;
     528             :                             }
     529             :                         } else {
     530             :                             // Shading surface is nearest hit
     531         484 :                             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPtSolRefl(RayNum, RecPtNum) =
     532         484 :                                 state.dataSurface->SurfShadowDiffuseSolRefl(NearestHitSurfNum);
     533             :                         }
     534             :                     } else {
     535             :                         // No obstructions were hit by this ray
     536       16096 :                         state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPtSurfNum(RayNum, RecPtNum) = 0;
     537             :                         // If ray is going downward find the hit point on the ground plane if the receiving point
     538             :                         // is above ground level; note that GroundLevelZ is <= 0.0
     539       23346 :                         if (RayVec(3) < 0.0 &&
     540        7250 :                             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPt(RecPtNum).z > state.dataSurface->GroundLevelZ) {
     541             :                             // Ray hits ground
     542        7250 :                             Alfa = std::acos(-RayVec.z);
     543        7250 :                             Beta = std::atan2(RayVec.y, RayVec.x);
     544        7250 :                             HorDis = (RecPt.z - state.dataSurface->GroundLevelZ) * std::tan(Alfa);
     545        7250 :                             GroundHitPt.z = state.dataSurface->GroundLevelZ;
     546        7250 :                             GroundHitPt.x = RecPt.x + HorDis * std::cos(Beta);
     547        7250 :                             GroundHitPt.y = RecPt.y + HorDis * std::sin(Beta);
     548        7250 :                             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPt(RayNum, RecPtNum) = GroundHitPt;
     549        7250 :                             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPtSurfNum(RayNum, RecPtNum) = -1;
     550        7250 :                             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPtHitPtDis(RayNum, RecPtNum) =
     551        7250 :                                 (RecPt(3) - state.dataSurface->GroundLevelZ) / (-RayVec(3));
     552        7250 :                             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPtSolRefl(RayNum, RecPtNum) =
     553        7250 :                                 state.dataEnvrn->GndReflectance;
     554        7250 :                             state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).HitPtNormVec(RayNum, RecPtNum) = unit_z;
     555             :                         } // End of check if ray hits ground
     556             :                     }     // End of check if obstruction hit
     557             :                 }         // End of RayNum loop
     558             :             }             // End of receiving point loop
     559             :         }                 // End of receiving surface loop
     560             :     }
     561             : 
     562             :     //=====================================================================================================
     563             : 
     564          54 :     void CalcBeamSolDiffuseReflFactors(EnergyPlusData &state)
     565             :     {
     566             : 
     567             :         // SUBROUTINE INFORMATION:
     568             :         //       AUTHOR         Fred Winkelmann
     569             :         //       DATE WRITTEN   September 2003
     570             :         //       MODIFIED       TH 4/6/2010, fixed CR 7872
     571             :         //       RE-ENGINEERED  B. Griffith, October 2012, for timestep integrated solar.
     572             : 
     573             :         // PURPOSE OF THIS SUBROUTINE:
     574             :         // manage calculations for factors for irradiance on exterior heat transfer surfaces due to
     575             :         // beam-to-diffuse solar reflection from obstructions and ground.
     576             : 
     577             :         // METHODOLOGY EMPLOYED: call worker routine depending on solar calculation method
     578             : 
     579          54 :         if (!state.dataSysVars->DetailedSolarTimestepIntegration) {
     580          54 :             if (state.dataGlobal->BeginSimFlag) {
     581           9 :                 DisplayString(state, "Calculating Beam-to-Diffuse Exterior Solar Reflection Factors");
     582             :             } else {
     583          45 :                 DisplayString(state, "Updating Beam-to-Diffuse Exterior Solar Reflection Factors");
     584             :             }
     585          54 :             state.dataSurface->SurfReflFacBmToDiffSolObs = 0.0;
     586          54 :             state.dataSurface->SurfReflFacBmToDiffSolGnd = 0.0;
     587        1350 :             for (state.dataSolarReflectionManager->IHr = 1; state.dataSolarReflectionManager->IHr <= 24; ++state.dataSolarReflectionManager->IHr) {
     588        1296 :                 FigureBeamSolDiffuseReflFactors(state, state.dataSolarReflectionManager->IHr);
     589             :             }    // End of IHr loop
     590             :         } else { // timestep integrated solar, use current hour of day
     591           0 :             state.dataSurface->SurfReflFacBmToDiffSolObs(state.dataGlobal->HourOfDay, {1, state.dataSurface->TotSurfaces}) = 0.0;
     592           0 :             state.dataSurface->SurfReflFacBmToDiffSolGnd(state.dataGlobal->HourOfDay, {1, state.dataSurface->TotSurfaces}) = 0.0;
     593           0 :             FigureBeamSolDiffuseReflFactors(state, state.dataGlobal->HourOfDay);
     594             :         }
     595          54 :     }
     596             : 
     597        1296 :     void FigureBeamSolDiffuseReflFactors(EnergyPlusData &state, int const iHour)
     598             :     {
     599             : 
     600             :         // SUBROUTINE INFORMATION:
     601             :         //       AUTHOR         Fred Winkelmann, derived from original CalcBeamSolDiffuseReflFactors
     602             :         //       DATE WRITTEN   September 2003
     603             :         //       MODIFIED       na
     604             :         //       RE-ENGINEERED  B. Griffith, October 2012, revised for timestep integrated solar
     605             : 
     606             :         // PURPOSE OF THIS SUBROUTINE:
     607             :         // Calculates factors for irradiance on exterior heat transfer surfaces due to
     608             :         // beam-to-diffuse solar reflection from obstructions and ground.
     609             : 
     610        2592 :         Array1D<Real64> ReflBmToDiffSolObs(state.dataSurface->MaxRecPts); // Irradiance at a receiving point for
     611             :         // beam solar diffusely reflected from obstructions, divided by
     612             :         // beam normal irradiance
     613        2592 :         Array1D<Real64> ReflBmToDiffSolGnd(state.dataSurface->MaxRecPts); // Irradiance at a receiving point for
     614             :         // beam solar diffusely reflected from the ground, divided by
     615             :         // beam normal irradiance
     616             :         bool hit; // True iff obstruction is hit
     617        1296 :         ReflBmToDiffSolObs = 0.0;
     618        1296 :         ReflBmToDiffSolGnd = 0.0;
     619             : 
     620             :         // Unit vector to sun
     621        1296 :         state.dataSolarReflectionManager->SunVec = state.dataSurface->SurfSunCosHourly(iHour);
     622             : 
     623             :         // loop through each surface that can receive beam solar reflected as diffuse solar from other surfaces
     624        7392 :         for (state.dataSolarReflectionManager->RecSurfNum = 1;
     625        7392 :              state.dataSolarReflectionManager->RecSurfNum <= state.dataSolarReflectionManager->TotSolReflRecSurf;
     626        6096 :              ++state.dataSolarReflectionManager->RecSurfNum) {
     627        6096 :             state.dataSolarReflectionManager->SurfNum =
     628        6096 :                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum).SurfNum;
     629             : 
     630       30480 :             for (state.dataSolarReflectionManager->RecPtNum = 1;
     631       60960 :                  state.dataSolarReflectionManager->RecPtNum <=
     632       30480 :                  state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum).NumRecPts;
     633       24384 :                  ++state.dataSolarReflectionManager->RecPtNum) {
     634       24384 :                 ReflBmToDiffSolObs(state.dataSolarReflectionManager->RecPtNum) = 0.0;
     635       24384 :                 ReflBmToDiffSolGnd(state.dataSolarReflectionManager->RecPtNum) = 0.0;
     636             : 
     637     2218944 :                 for (state.dataSolarReflectionManager->RayNum = 1;
     638     4437888 :                      state.dataSolarReflectionManager->RayNum <=
     639     2218944 :                      state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum).NumReflRays;
     640     2194560 :                      ++state.dataSolarReflectionManager->RayNum) {
     641     2194560 :                     state.dataSolarReflectionManager->HitPtSurfNum =
     642     2194560 :                         state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum)
     643     2194560 :                             .HitPtSurfNum(state.dataSolarReflectionManager->RayNum, state.dataSolarReflectionManager->RecPtNum);
     644             : 
     645             :                     // Skip rays that do not hit an obstruction or ground.
     646             :                     // (Note that if a downgoing ray does not hit an obstruction it will have HitPtSurfNum = 0
     647             :                     // if the receiving point is below ground level (see subr. InitSolReflRecSurf); this means
     648             :                     // that a below-ground-level receiving point receives no ground-reflected radiation although
     649             :                     // it is allowed to receive obstruction-reflected solar radiation and direct (unreflected)
     650             :                     // beam and sky solar radiation. As far as reflected solar is concerned, the program does
     651             :                     // not handle a sloped ground plane or a horizontal ground plane whose level is different
     652             :                     // from one side of the building to another.)
     653     2194560 :                     if (state.dataSolarReflectionManager->HitPtSurfNum == 0)
     654     1144992 :                         continue; // Ray hits sky or obstruction with receiving pt. below ground level
     655             : 
     656     1049568 :                     if (state.dataSolarReflectionManager->HitPtSurfNum > 0) {
     657             :                         // Skip rays that hit a daylighting shelf, from which solar reflection is calculated separately.
     658       91776 :                         if (state.dataSurface->SurfDaylightingShelfInd(state.dataSolarReflectionManager->HitPtSurfNum) > 0) continue;
     659             : 
     660             :                         // Skip rays that hit a window
     661             :                         // If hit point's surface is a window or glass door go to next ray since it is assumed for now
     662             :                         // that windows have only beam-to-beam, not beam-to-diffuse, reflection
     663             :                         // TH 3/29/2010. Code modified and moved
     664      171072 :                         if (state.dataSurface->Surface(state.dataSolarReflectionManager->HitPtSurfNum).Class == SurfaceClass::Window ||
     665       79296 :                             state.dataSurface->Surface(state.dataSolarReflectionManager->HitPtSurfNum).Class == SurfaceClass::GlassDoor)
     666       12480 :                             continue;
     667             : 
     668             :                         // Skip rays that hit non-sunlit surface. Assume first time step of the hour.
     669       79296 :                         state.dataSolarReflectionManager->SunLitFract =
     670       79296 :                             state.dataHeatBal->SurfSunlitFrac(iHour, 1, state.dataSolarReflectionManager->HitPtSurfNum);
     671             : 
     672             :                         // If hit point's surface is not sunlit go to next ray
     673             :                         // TH 3/25/2010. why limit to HeatTransSurf? shading surfaces should also apply
     674             :                         // IF(Surface(HitPtSurfNum)%HeatTransSurf .AND. SunLitFract < 0.01d0) CYCLE
     675       79296 :                         if (state.dataSolarReflectionManager->SunLitFract < 0.01) continue;
     676             : 
     677             :                         // TH 3/26/2010. If the hit point falls into the shadow even though SunLitFract > 0, can Cycle.
     678             :                         //  This cannot be done now, therefore there are follow-up checks of blocking sun ray
     679             :                         //   from the hit point.
     680             : 
     681             :                         // TH 3/29/2010. Code modified and moved up
     682             :                         // If hit point's surface is a window go to next ray since it is assumed for now
     683             :                         // that windows have only beam-to-beam, not beam-to-diffuse, reflection
     684             :                         // IF(Surface(HitPtSurfNum)%Construction > 0) THEN
     685             :                         //  IF(Construct(Surface(HitPtSurfNum)%Construction)%TypeIsWindow) CYCLE
     686             :                         // END IF
     687             :                     }
     688             : 
     689             :                     // Does an obstruction block the vector from this ray's hit point to the sun?
     690     1928220 :                     state.dataSolarReflectionManager->OriginThisRay =
     691      964110 :                         state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum)
     692     1928220 :                             .HitPt(state.dataSolarReflectionManager->RayNum, state.dataSolarReflectionManager->RecPtNum);
     693             : 
     694             :                     // Note: if sun is in back of hit surface relative to receiving point, CosIncBmAtHitPt will be < 0
     695      964110 :                     state.dataSolarReflectionManager->CosIncBmAtHitPt =
     696     1928220 :                         dot(state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum)
     697      964110 :                                 .HitPtNormVec(state.dataSolarReflectionManager->RayNum, state.dataSolarReflectionManager->RecPtNum),
     698      964110 :                             state.dataSolarReflectionManager->SunVec);
     699      964110 :                     if (state.dataSolarReflectionManager->CosIncBmAtHitPt <= 0.0) continue;
     700             : 
     701             :                     // CR 7872 - TH 4/6/2010. The shading surfaces should point to the receiveing heat transfer surface
     702             :                     //  according to the the right hand rule. If user inputs do not follow the rule, use the following
     703             :                     //  code to check the mirrored shading surface
     704      226753 :                     if (state.dataSolarReflectionManager->HitPtSurfNum > 0) {
     705        5593 :                         if (state.dataSurface->Surface(state.dataSolarReflectionManager->HitPtSurfNum).IsShadowing) {
     706        3328 :                             if (state.dataSolarReflectionManager->HitPtSurfNum + 1 < state.dataSurface->TotSurfaces) {
     707        3328 :                                 if (state.dataSurface->Surface(state.dataSolarReflectionManager->HitPtSurfNum + 1).IsShadowing &&
     708           0 :                                     state.dataSurface->Surface(state.dataSolarReflectionManager->HitPtSurfNum + 1).MirroredSurf) {
     709             :                                     // Check whether the sun is behind the mirrored shading surface
     710           0 :                                     state.dataSolarReflectionManager->CosIncBmAtHitPt2 =
     711           0 :                                         dot(state.dataSurface->Surface(state.dataSolarReflectionManager->HitPtSurfNum + 1).OutNormVec,
     712           0 :                                             state.dataSolarReflectionManager->SunVec);
     713           0 :                                     if (state.dataSolarReflectionManager->CosIncBmAtHitPt2 >= 0.0) continue;
     714             :                                 }
     715             :                             }
     716             :                         }
     717             :                     }
     718             : 
     719             :                     // TH 3/25/2010. CR 7872. Seems should loop over all possible obstructions for the HitPtSurfNum
     720             :                     //  rather than RecSurfNum, because if the HitPtSurfNum is a shading surface,
     721             :                     //  it does not belong to SolReflRecSurf which only contain heat transfer surfaces
     722             :                     //  that can receive reflected solar (ExtSolar = True)!
     723             : 
     724             :                     // To speed up, ideally should store all possible shading surfaces for the HitPtSurfNum
     725             :                     //  obstruction surface in the SolReflSurf(HitPtSurfNum)%PossibleObsSurfNums(loop) array as well
     726      226753 :                     hit = false;
     727     2455219 :                     for (state.dataSolarReflectionManager->ObsSurfNum = 1;
     728     2455219 :                          state.dataSolarReflectionManager->ObsSurfNum <= state.dataSurface->TotSurfaces;
     729     2228466 :                          ++state.dataSolarReflectionManager->ObsSurfNum) {
     730             :                         //        DO loop = 1,SolReflRecSurf(RecSurfNum)%NumPossibleObs
     731             :                         //          ObsSurfNum = SolReflRecSurf(RecSurfNum)%PossibleObsSurfNums(loop)
     732             : 
     733             :                         // CR 8959 -- The other side of a mirrored surface cannot obstruct the mirrored surface
     734     2281532 :                         if (state.dataSolarReflectionManager->HitPtSurfNum > 0) {
     735       99812 :                             if (state.dataSurface->Surface(state.dataSolarReflectionManager->HitPtSurfNum).MirroredSurf) {
     736       79872 :                                 if (state.dataSolarReflectionManager->ObsSurfNum == state.dataSolarReflectionManager->HitPtSurfNum - 1) continue;
     737             :                             }
     738             :                         }
     739             : 
     740             :                         // skip the hit surface
     741     2278204 :                         if (state.dataSolarReflectionManager->ObsSurfNum == state.dataSolarReflectionManager->HitPtSurfNum) continue;
     742             : 
     743             :                         // skip mirrored surfaces
     744     2274072 :                         if (state.dataSurface->Surface(state.dataSolarReflectionManager->ObsSurfNum).MirroredSurf) continue;
     745             :                         // IF(Surface(ObsSurfNum)%ShadowingSurf .AND. Surface(ObsSurfNum)%Name(1:3) == 'Mir') THEN
     746             :                         //  CYCLE
     747             :                         // ENDIF
     748             : 
     749             :                         // skip interior surfaces
     750     2174495 :                         if (state.dataSurface->Surface(state.dataSolarReflectionManager->ObsSurfNum).ExtBoundCond >= 1) continue;
     751             : 
     752             :                         // For now it is assumed that obstructions that are shading surfaces are opaque.
     753             :                         // An improvement here would be to allow these to have transmittance.
     754     5182404 :                         PierceSurface(state,
     755     1295601 :                                       state.dataSolarReflectionManager->ObsSurfNum,
     756     1295601 :                                       state.dataSolarReflectionManager->OriginThisRay,
     757     1295601 :                                       state.dataSolarReflectionManager->SunVec,
     758     1295601 :                                       state.dataSolarReflectionManager->ObsHitPt,
     759             :                                       hit);
     760     1295601 :                         if (hit) break; // An obstruction was hit
     761             :                     }
     762      226753 :                     if (hit) continue; // Sun does not reach this ray's hit point
     763             : 
     764             :                     // Sun reaches this ray's hit point; get beam-reflected diffuse radiance at hit point for
     765             :                     // unit beam normal solar
     766             : 
     767             :                     // CosIncBmAtHitPt = DOT_PRODUCT(SolReflRecSurf(RecSurfNum)%HitPtNormVec(RecPtNum,RayNum),SunVec)
     768             :                     // Note: if sun is in back of hit surface relative to receiving point, CosIncBmAtHitPt will be < 0
     769             :                     // and use of MAX in following gives zero beam solar reflecting at hit point.
     770             :                     // BmReflSolRadiance = MAX(0.0d0,CosIncBmAtHitPt)*SolReflRecSurf(RecSurfNum)%HitPtSolRefl(RecPtNum,RayNum)
     771             : 
     772      173687 :                     state.dataSolarReflectionManager->BmReflSolRadiance =
     773      347374 :                         state.dataSolarReflectionManager->CosIncBmAtHitPt *
     774      173687 :                         state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum)
     775      173687 :                             .HitPtSolRefl(state.dataSolarReflectionManager->RayNum, state.dataSolarReflectionManager->RecPtNum);
     776             : 
     777      173687 :                     if (state.dataSolarReflectionManager->BmReflSolRadiance > 0.0) {
     778             :                         // Contribution to reflection factor from this hit point
     779      173687 :                         if (state.dataSolarReflectionManager->HitPtSurfNum > 0) {
     780             :                             // Ray hits an obstruction
     781        4079 :                             state.dataSolarReflectionManager->dReflBeamToDiffSol =
     782        8158 :                                 state.dataSolarReflectionManager->BmReflSolRadiance *
     783        4079 :                                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum)
     784        8158 :                                     .dOmegaRay(state.dataSolarReflectionManager->RayNum) *
     785        4079 :                                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum)
     786        8158 :                                     .CosIncAngRay(state.dataSolarReflectionManager->RayNum) /
     787             :                                 DataGlobalConstants::Pi;
     788        4079 :                             ReflBmToDiffSolObs(state.dataSolarReflectionManager->RecPtNum) += state.dataSolarReflectionManager->dReflBeamToDiffSol;
     789             :                         } else {
     790             :                             // Ray hits ground (in this case we do not multiply by BmReflSolRadiance since
     791             :                             // ground reflectance and cos of incidence angle of sun on
     792             :                             // ground is taken into account later when SurfReflFacBmToDiffSolGnd is used)
     793      169608 :                             state.dataSolarReflectionManager->dReflBeamToDiffSol =
     794      169608 :                                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum)
     795      339216 :                                     .dOmegaRay(state.dataSolarReflectionManager->RayNum) *
     796      169608 :                                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum)
     797      339216 :                                     .CosIncAngRay(state.dataSolarReflectionManager->RayNum) /
     798             :                                 DataGlobalConstants::Pi;
     799      169608 :                             ReflBmToDiffSolGnd(state.dataSolarReflectionManager->RecPtNum) += state.dataSolarReflectionManager->dReflBeamToDiffSol;
     800             :                         }
     801             :                     }
     802             :                 } // End of loop over rays from receiving point
     803             :             }     // End of loop over receiving points
     804             : 
     805             :             // Average over receiving points
     806        6096 :             state.dataSurface->SurfReflFacBmToDiffSolObs(iHour, state.dataSolarReflectionManager->SurfNum) = 0.0;
     807        6096 :             state.dataSurface->SurfReflFacBmToDiffSolGnd(iHour, state.dataSolarReflectionManager->SurfNum) = 0.0;
     808        6096 :             state.dataSolarReflectionManager->NumRecPts =
     809        6096 :                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->RecSurfNum).NumRecPts;
     810       30480 :             for (state.dataSolarReflectionManager->RecPtNum = 1;
     811       30480 :                  state.dataSolarReflectionManager->RecPtNum <= state.dataSolarReflectionManager->NumRecPts;
     812       24384 :                  ++state.dataSolarReflectionManager->RecPtNum) {
     813       24384 :                 state.dataSurface->SurfReflFacBmToDiffSolObs(iHour, state.dataSolarReflectionManager->SurfNum) +=
     814       24384 :                     ReflBmToDiffSolObs(state.dataSolarReflectionManager->RecPtNum);
     815       24384 :                 state.dataSurface->SurfReflFacBmToDiffSolGnd(iHour, state.dataSolarReflectionManager->SurfNum) +=
     816       24384 :                     ReflBmToDiffSolGnd(state.dataSolarReflectionManager->RecPtNum);
     817             :             }
     818        6096 :             state.dataSurface->SurfReflFacBmToDiffSolObs(iHour, state.dataSolarReflectionManager->SurfNum) /=
     819        6096 :                 state.dataSolarReflectionManager->NumRecPts;
     820        6096 :             state.dataSurface->SurfReflFacBmToDiffSolGnd(iHour, state.dataSolarReflectionManager->SurfNum) /=
     821        6096 :                 state.dataSolarReflectionManager->NumRecPts;
     822             : 
     823             :             // Do not allow SurfReflFacBmToDiffSolGnd to exceed the surface's unobstructed ground view factor
     824        6096 :             state.dataSurface->SurfReflFacBmToDiffSolGnd(iHour, state.dataSolarReflectionManager->SurfNum) =
     825        6096 :                 min(0.5 * (1.0 - state.dataSurface->Surface(state.dataSolarReflectionManager->SurfNum).CosTilt),
     826        6096 :                     state.dataSurface->SurfReflFacBmToDiffSolGnd(iHour, state.dataSolarReflectionManager->SurfNum));
     827             :             // Note: the above factors are dimensionless; they are equal to
     828             :             // (W/m2 reflected solar incident on SurfNum)/(W/m2 beam normal solar)
     829             :         } // End of loop over receiving surfaces
     830        1296 :     }
     831             : 
     832             :     //=================================================================================================
     833             : 
     834          54 :     void CalcBeamSolSpecularReflFactors(EnergyPlusData &state)
     835             :     {
     836             : 
     837             :         // SUBROUTINE INFORMATION:
     838             :         //       AUTHOR         Fred Winkelmann
     839             :         //       DATE WRITTEN   September 2003
     840             :         //       MODIFIED       na
     841             :         //       RE-ENGINEERED  B. Griffith, October 2012, for timestep integrated solar
     842             : 
     843             :         // PURPOSE OF THIS SUBROUTINE:
     844             :         // Manage calculation of factors for beam solar irradiance on exterior heat transfer surfaces due to
     845             :         // specular (beam-to-beam) reflection from obstructions such as a highly-glazed neighboring
     846             :         // building.
     847             : 
     848             :         // METHODOLOGY EMPLOYED:
     849             :         // call worker routine as appropriate
     850             : 
     851          54 :         if (!state.dataSysVars->DetailedSolarTimestepIntegration) {
     852          54 :             if (state.dataGlobal->BeginSimFlag) {
     853           9 :                 DisplayString(state, "Calculating Beam-to-Beam Exterior Solar Reflection Factors");
     854             :             } else {
     855          45 :                 DisplayString(state, "Updating Beam-to-Beam Exterior Solar Reflection Factors");
     856             :             }
     857          54 :             state.dataSurface->SurfReflFacBmToBmSolObs = 0.0;
     858          54 :             state.dataSurface->SurfCosIncAveBmToBmSolObs = 0.0;
     859        1350 :             for (state.dataSolarReflectionManager->NumHr = 1; state.dataSolarReflectionManager->NumHr <= 24;
     860        1296 :                  ++state.dataSolarReflectionManager->NumHr) {
     861        1296 :                 FigureBeamSolSpecularReflFactors(state, state.dataSolarReflectionManager->NumHr);
     862             :             }    // End of NumHr loop
     863             :         } else { // timestep integrated solar, use current hour of day
     864           0 :             state.dataSurface->SurfReflFacBmToBmSolObs(state.dataGlobal->HourOfDay, {1, state.dataSurface->TotSurfaces}) = 0.0;
     865           0 :             state.dataSurface->SurfCosIncAveBmToBmSolObs(state.dataGlobal->HourOfDay, {1, state.dataSurface->TotSurfaces}) = 0.0;
     866           0 :             FigureBeamSolSpecularReflFactors(state, state.dataGlobal->HourOfDay);
     867             :         }
     868          54 :     }
     869             : 
     870        1296 :     void FigureBeamSolSpecularReflFactors(EnergyPlusData &state, int const iHour)
     871             :     {
     872             : 
     873             :         // SUBROUTINE INFORMATION:
     874             :         //       AUTHOR         Fred Winkelmann
     875             :         //       DATE WRITTEN   September 2003
     876             :         //       MODIFIED       na
     877             :         //       RE-ENGINEERED  B. Griffith, October 2012, for timestep integrated solar
     878             : 
     879             :         // PURPOSE OF THIS SUBROUTINE:
     880             :         // Calculates factors for beam solar irradiance on exterior heat transfer surfaces due to
     881             :         // specular (beam-to-beam) reflection from obstructions such as a highly-glazed neighboring
     882             :         // building. Specular reflection can occur from shading surfaces with non-zero specular
     883             :         // reflectance and from exterior windows of the building (in calculating reflection from
     884             :         // these windows, they are assumed to have no shades or blinds).
     885             :         // Reflection from the ground and opaque building surfaces is assumed to be totally diffuse,
     886             :         // i.e. these surfaces has no specular reflection component.
     887             : 
     888             :         // METHODOLOGY EMPLOYED:
     889             :         // <description>
     890             : 
     891             :         // REFERENCES:
     892             :         // na
     893             : 
     894             :         // Using/Aliasing
     895             :         using General::POLYF;
     896             : 
     897             :         // Locals
     898             :         // SUBROUTINE LOCAL VARIABLE DECLARATIONS:
     899        1648 :         Array1D<Real64> ReflBmToDiffSolObs(state.dataSurface->MaxRecPts); // Irradiance at a receiving point for
     900             :         // beam solar diffusely reflected from obstructions, divided by
     901             :         // beam normal irradiance
     902             :         // unused  INTEGER           :: RayNum               =0   ! Ray number
     903             :         bool hitRefl;                                                   // True iff reflecting surface is hit
     904             :         bool hitObs;                                                    // True iff obstruction is hit
     905             :         bool hitObsRefl;                                                // True iff obstruction hit between rec. pt. and reflection point
     906        1648 :         Array1D<Real64> ReflBmToBmSolObs(state.dataSurface->MaxRecPts); // Irradiance at a receiving point for
     907             :         // beam solar specularly reflected from obstructions, divided by
     908             :         // beam normal irradiance
     909             :         Real64 ReflDistanceSq; // Distance squared from receiving point to hit point on a reflecting surface (m)
     910             :         Real64 ReflDistance;   // Distance from receiving point to hit point on a reflecting surface (m)
     911        1648 :         Array1D<Real64> ReflFacTimesCosIncSum(state.dataSurface->MaxRecPts); // Sum of ReflFac times CosIncAngRefl
     912             : 
     913        1296 :         ReflBmToDiffSolObs = 0.0;
     914        1296 :         ReflFacTimesCosIncSum = 0.0;
     915             : 
     916        1296 :         if (state.dataSurface->SurfSunCosHourly(iHour)(3) < DataEnvironment::SunIsUpValue) return; // Skip if sun is below horizon
     917             : 
     918             :         // Unit vector to sun
     919         352 :         state.dataSolarReflectionManager->SunVect = state.dataSurface->SurfSunCosHourly(iHour);
     920             : 
     921        1740 :         for (int RecSurfNum = 1; RecSurfNum <= state.dataSolarReflectionManager->TotSolReflRecSurf; ++RecSurfNum) {
     922             :             int const SurfNum =
     923        1388 :                 state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).SurfNum; // Heat transfer surface number corresponding to RecSurfNum
     924        1388 :             if (state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumPossibleObs > 0) {
     925         384 :                 ReflBmToBmSolObs = 0.0;
     926         384 :                 ReflFacTimesCosIncSum = 0.0;
     927         384 :                 int const NumRecPts = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumRecPts;
     928             :                 // Find possible reflecting surfaces for this receiving surface
     929        1200 :                 for (int loop = 1, loop_end = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumPossibleObs; loop <= loop_end; ++loop) {
     930             :                     int const ReflSurfNum =
     931         816 :                         state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).PossibleObsSurfNums(loop); // Reflecting surface number
     932             :                     // Keep windows; keep shading surfaces with specular reflectance
     933        2112 :                     if ((state.dataSurface->Surface(ReflSurfNum).Class == SurfaceClass::Window && state.dataSurface->Surface(ReflSurfNum).ExtSolar) ||
     934        1248 :                         (state.dataSurface->SurfShadowGlazingFrac(ReflSurfNum) > 0.0 && state.dataSurface->Surface(ReflSurfNum).IsShadowing)) {
     935             :                         // Skip if window and not sunlit
     936         576 :                         if (state.dataSurface->Surface(ReflSurfNum).Class == SurfaceClass::Window &&
     937          48 :                             state.dataHeatBal->SurfSunlitFrac(iHour, 1, ReflSurfNum) < 0.01)
     938          26 :                             continue;
     939             :                         // Check if sun is in front of this reflecting surface.
     940         502 :                         state.dataSolarReflectionManager->ReflNorm = state.dataSurface->Surface(ReflSurfNum).OutNormVec;
     941         502 :                         state.dataSolarReflectionManager->CosIncAngRefl =
     942         502 :                             dot(state.dataSolarReflectionManager->SunVect, state.dataSolarReflectionManager->ReflNorm);
     943         502 :                         if (state.dataSolarReflectionManager->CosIncAngRefl < 0.0) continue;
     944             : 
     945             :                         // Get sun position unit vector for mirror image of sun in reflecting surface
     946         522 :                         state.dataSolarReflectionManager->SunVecMir =
     947         522 :                             state.dataSolarReflectionManager->SunVect -
     948         522 :                             2.0 * dot(state.dataSolarReflectionManager->SunVect, state.dataSolarReflectionManager->ReflNorm) *
     949         522 :                                 state.dataSolarReflectionManager->ReflNorm;
     950             :                         // Angle of incidence of reflected beam on receiving surface
     951         261 :                         state.dataSolarReflectionManager->CosIncAngRec =
     952         261 :                             dot(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec, state.dataSolarReflectionManager->SunVecMir);
     953         261 :                         if (state.dataSolarReflectionManager->CosIncAngRec <= 0.0) continue;
     954         870 :                         for (int RecPtNum = 1; RecPtNum <= NumRecPts; ++RecPtNum) {
     955             :                             // See if ray from receiving point to mirrored sun hits the reflecting surface
     956         696 :                             state.dataSolarReflectionManager->RecPt = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).RecPt(RecPtNum);
     957        2088 :                             PierceSurface(state,
     958             :                                           ReflSurfNum,
     959         696 :                                           state.dataSolarReflectionManager->RecPt,
     960         696 :                                           state.dataSolarReflectionManager->SunVecMir,
     961         696 :                                           state.dataSolarReflectionManager->HitPtRefl,
     962             :                                           hitRefl);
     963         696 :                             if (hitRefl) { // Reflecting surface was hit
     964         122 :                                 ReflDistanceSq =
     965         122 :                                     distance_squared(state.dataSolarReflectionManager->HitPtRefl, state.dataSolarReflectionManager->RecPt);
     966         122 :                                 ReflDistance = std::sqrt(ReflDistanceSq);
     967             :                                 // Determine if ray from receiving point to hit point is obstructed
     968         122 :                                 hitObsRefl = false;
     969         383 :                                 for (int loop2 = 1, loop2_end = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NumPossibleObs;
     970         383 :                                      loop2 <= loop2_end;
     971             :                                      ++loop2) {
     972         266 :                                     int const ObsSurfNum = state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).PossibleObsSurfNums(loop2);
     973         266 :                                     if (ObsSurfNum == ReflSurfNum || ObsSurfNum == state.dataSurface->Surface(ReflSurfNum).BaseSurf) continue;
     974         432 :                                     PierceSurface(state,
     975             :                                                   ObsSurfNum,
     976         144 :                                                   state.dataSolarReflectionManager->RecPt,
     977         144 :                                                   state.dataSolarReflectionManager->SunVecMir,
     978             :                                                   ReflDistance,
     979         144 :                                                   state.dataSolarReflectionManager->HitPtObs,
     980             :                                                   hitObs); // ReflDistance cutoff added
     981         144 :                                     if (hitObs) {          // => Could skip distance check (unless < vs <= ReflDistance really matters)
     982          84 :                                         if (distance_squared(state.dataSolarReflectionManager->HitPtObs, state.dataSolarReflectionManager->RecPt) <
     983             :                                             ReflDistanceSq) {
     984           5 :                                             hitObsRefl = true;
     985           5 :                                             break;
     986             :                                         }
     987             :                                     }
     988             :                                 }
     989         122 :                                 if (hitObsRefl) continue; // Obstruction closer than reflection pt. was hit; go to next rec. pt.
     990             :                                 // There is no obstruction for this ray between rec. pt. and hit point on reflecting surface.
     991             :                                 // See if ray from hit pt. on reflecting surface to original (unmirrored) sun position is obstructed
     992         117 :                                 hitObs = false;
     993         117 :                                 if (state.dataSurface->Surface(ReflSurfNum).Class == SurfaceClass::Window) { // Reflecting surface is a window
     994             :                                     // Receiving surface number for this window
     995           0 :                                     int const ReflSurfRecNum = state.dataSurface->SurfShadowRecSurfNum(
     996           0 :                                         ReflSurfNum); // Receiving surface number corresponding to a reflecting surface number
     997           0 :                                     if (ReflSurfRecNum > 0) {
     998             :                                         // Loop over possible obstructions for this window
     999           0 :                                         for (int loop2 = 1,
    1000           0 :                                                  loop2_end = state.dataSolarReflectionManager->SolReflRecSurf(ReflSurfRecNum).NumPossibleObs;
    1001           0 :                                              loop2 <= loop2_end;
    1002             :                                              ++loop2) {
    1003             :                                             int const ObsSurfNum =
    1004           0 :                                                 state.dataSolarReflectionManager->SolReflRecSurf(ReflSurfRecNum).PossibleObsSurfNums(loop2);
    1005           0 :                                             PierceSurface(state,
    1006             :                                                           ObsSurfNum,
    1007           0 :                                                           state.dataSolarReflectionManager->HitPtRefl,
    1008           0 :                                                           state.dataSolarReflectionManager->SunVect,
    1009           0 :                                                           state.dataSolarReflectionManager->HitPtObs,
    1010             :                                                           hitObs);
    1011           0 :                                             if (hitObs) break;
    1012             :                                         }
    1013             :                                     }
    1014             :                                 } else { // Reflecting surface is a building shade
    1015        1521 :                                     for (int ObsSurfNum : state.dataSurface->AllShadowPossObstrSurfaceList) {
    1016        1404 :                                         if (ObsSurfNum == ReflSurfNum) continue;
    1017             : 
    1018             :                                         // TH2 CR8959 -- Skip mirrored surfaces
    1019        1404 :                                         if (state.dataSurface->Surface(ObsSurfNum).MirroredSurf) continue;
    1020             :                                         // TH2 CR8959 -- The other side of a mirrored surface cannot obstruct the mirrored surface
    1021        1404 :                                         if (state.dataSurface->Surface(ReflSurfNum).MirroredSurf) {
    1022        1404 :                                             if (ObsSurfNum == ReflSurfNum - 1) continue;
    1023             :                                         }
    1024             : 
    1025        3861 :                                         PierceSurface(state,
    1026             :                                                       ObsSurfNum,
    1027        1287 :                                                       state.dataSolarReflectionManager->HitPtRefl,
    1028        1287 :                                                       state.dataSolarReflectionManager->SunVect,
    1029        1287 :                                                       state.dataSolarReflectionManager->HitPtObs,
    1030             :                                                       hitObs);
    1031        1287 :                                         if (hitObs) break;
    1032             :                                     }
    1033             :                                 }
    1034             : 
    1035         117 :                                 if (hitObs) continue; // Obstruction hit between reflection hit point and sun; go to next receiving pt.
    1036             : 
    1037             :                                 // No obstructions. Calculate reflected beam irradiance at receiving pt. from this reflecting surface.
    1038         117 :                                 state.dataSolarReflectionManager->SpecReflectance = 0.0;
    1039         117 :                                 if (state.dataSurface->Surface(ReflSurfNum).Class == SurfaceClass::Window) {
    1040           0 :                                     state.dataSolarReflectionManager->ConstrNumRefl = state.dataSurface->Surface(ReflSurfNum).Construction;
    1041           0 :                                     state.dataSolarReflectionManager->SpecReflectance = POLYF(
    1042           0 :                                         std::abs(state.dataSolarReflectionManager->CosIncAngRefl),
    1043           0 :                                         state.dataConstruction->Construct(state.dataSolarReflectionManager->ConstrNumRefl).ReflSolBeamFrontCoef);
    1044             :                                 }
    1045         234 :                                 if (state.dataSurface->Surface(ReflSurfNum).IsShadowing &&
    1046         117 :                                     state.dataSurface->SurfShadowGlazingConstruct(ReflSurfNum) > 0) {
    1047         117 :                                     state.dataSolarReflectionManager->ConstrNumRefl = state.dataSurface->SurfShadowGlazingConstruct(ReflSurfNum);
    1048         117 :                                     state.dataSolarReflectionManager->SpecReflectance =
    1049         234 :                                         state.dataSurface->SurfShadowGlazingFrac(ReflSurfNum) *
    1050         234 :                                         POLYF(
    1051         117 :                                             std::abs(state.dataSolarReflectionManager->CosIncAngRefl),
    1052         117 :                                             state.dataConstruction->Construct(state.dataSolarReflectionManager->ConstrNumRefl).ReflSolBeamFrontCoef);
    1053             :                                 }
    1054             :                                 // Angle of incidence of reflected beam on receiving surface
    1055         117 :                                 state.dataSolarReflectionManager->CosIncAngRec =
    1056         117 :                                     dot(state.dataSolarReflectionManager->SolReflRecSurf(RecSurfNum).NormVec,
    1057         117 :                                         state.dataSolarReflectionManager->SunVecMir);
    1058         117 :                                 state.dataSolarReflectionManager->ReflFac =
    1059         117 :                                     state.dataSolarReflectionManager->SpecReflectance * state.dataSolarReflectionManager->CosIncAngRec;
    1060             :                                 // Contribution to specular reflection factor
    1061         117 :                                 ReflBmToBmSolObs(RecPtNum) += state.dataSolarReflectionManager->ReflFac;
    1062         117 :                                 ReflFacTimesCosIncSum(RecPtNum) +=
    1063         117 :                                     state.dataSolarReflectionManager->ReflFac * state.dataSolarReflectionManager->CosIncAngRec;
    1064             :                             } // End of check if reflecting surface was hit
    1065             :                         }     // End of loop over receiving points
    1066             :                     }         // End of check if valid reflecting surface
    1067             :                 }             // End of loop over obstructing surfaces
    1068             :                 // Average over receiving points
    1069             : 
    1070        1920 :                 for (int RecPtNum = 1; RecPtNum <= NumRecPts; ++RecPtNum) {
    1071        1536 :                     if (ReflBmToBmSolObs(RecPtNum) != 0.0) {
    1072         117 :                         state.dataSolarReflectionManager->CosIncWeighted = ReflFacTimesCosIncSum(RecPtNum) / ReflBmToBmSolObs(RecPtNum);
    1073             :                     } else {
    1074        1419 :                         state.dataSolarReflectionManager->CosIncWeighted = 0.0;
    1075             :                     }
    1076        1536 :                     state.dataSurface->SurfCosIncAveBmToBmSolObs(iHour, SurfNum) += state.dataSolarReflectionManager->CosIncWeighted;
    1077        1536 :                     state.dataSurface->SurfReflFacBmToBmSolObs(iHour, SurfNum) += ReflBmToBmSolObs(RecPtNum);
    1078             :                 }
    1079         384 :                 state.dataSurface->SurfReflFacBmToBmSolObs(iHour, SurfNum) /= double(NumRecPts);
    1080         384 :                 state.dataSurface->SurfCosIncAveBmToBmSolObs(iHour, SurfNum) /= double(NumRecPts);
    1081             :             } // End of check if number of possible obstructions > 0
    1082             :         }     // End of loop over receiving surfaces
    1083             :     }
    1084             : 
    1085             :     //=================================================================================================
    1086             : 
    1087           9 :     void CalcSkySolDiffuseReflFactors(EnergyPlusData &state)
    1088             :     {
    1089             : 
    1090             :         // SUBROUTINE INFORMATION:
    1091             :         //       AUTHOR         Fred Winkelmann
    1092             :         //       DATE WRITTEN   October 2003
    1093             :         //       MODIFIED       na
    1094             :         //       RE-ENGINEERED  na
    1095             : 
    1096             :         // PURPOSE OF THIS SUBROUTINE:
    1097             :         // Calculates factors for irradiance on exterior heat transfer surfaces due to
    1098             :         // reflection of sky diffuse solar radiation from obstructions and ground.
    1099             : 
    1100             :         // Using/Aliasing
    1101             : 
    1102             :         // SUBROUTINE LOCAL VARIABLE DECLARATIONS:
    1103          18 :         Array1D<Real64> ReflSkySolObs(state.dataSurface->MaxRecPts); // Irradiance at a receiving point for sky diffuse solar
    1104             :         // reflected from obstructions, divided by unobstructed
    1105             :         // sky diffuse horizontal irradiance
    1106          18 :         Array1D<Real64> ReflSkySolGnd(state.dataSurface->MaxRecPts); // Irradiance at a receiving point for sky diffuse solar
    1107             :         // reflected from ground, divided by unobstructed
    1108             :         // sky diffuse horizontal irradiance
    1109             :         bool hitObs; // True iff obstruction is hit
    1110             : 
    1111           9 :         Real64 const DPhi(DataGlobalConstants::PiOvr2 / (AltAngStepsForSolReflCalc / 2.0));      // Altitude angle and increment (radians)
    1112           9 :         Real64 const DTheta(2.0 * DataGlobalConstants::Pi / (2.0 * AzimAngStepsForSolReflCalc)); // Azimuth increment (radians)
    1113             : 
    1114             :         // Pre-compute these constants
    1115             :         // Initialize the 0 index with dummy value so the iterators line up below
    1116          18 :         std::vector<Real64> sin_Phi({-1});
    1117          18 :         std::vector<Real64> cos_Phi({-1});
    1118          54 :         for (int IPhi = 1; IPhi <= (AltAngStepsForSolReflCalc / 2); ++IPhi) {
    1119          45 :             Real64 Phi((IPhi - 0.5) * DPhi); // Altitude angle and increment (radians)
    1120          45 :             sin_Phi.push_back(std::sin(Phi));
    1121          45 :             cos_Phi.push_back(std::cos(Phi));
    1122             :         }
    1123             : 
    1124          18 :         std::vector<Real64> sin_Theta({-1});
    1125          18 :         std::vector<Real64> cos_Theta({-1});
    1126         171 :         for (int ITheta = 1; ITheta <= 2 * AzimAngStepsForSolReflCalc; ++ITheta) {
    1127         162 :             Real64 Theta = (ITheta - 0.5) * DTheta; // Azimuth angle (radians)
    1128         162 :             sin_Theta.push_back(std::sin(Theta));
    1129         162 :             cos_Theta.push_back(std::cos(Theta));
    1130             :         }
    1131             : 
    1132           9 :         DisplayString(state, "Calculating Sky Diffuse Exterior Solar Reflection Factors");
    1133           9 :         ReflSkySolObs = 0.0;
    1134           9 :         ReflSkySolGnd = 0.0;
    1135             : 
    1136          56 :         for (state.dataSolarReflectionManager->iRecSurfNum = 1;
    1137          56 :              state.dataSolarReflectionManager->iRecSurfNum <= state.dataSolarReflectionManager->TotSolReflRecSurf;
    1138          47 :              ++state.dataSolarReflectionManager->iRecSurfNum) {
    1139          47 :             state.dataSolarReflectionManager->iSurfNum =
    1140          47 :                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum).SurfNum;
    1141         235 :             for (state.dataSolarReflectionManager->iRecPtNum = 1;
    1142         470 :                  state.dataSolarReflectionManager->iRecPtNum <=
    1143         235 :                  state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum).NumRecPts;
    1144         188 :                  ++state.dataSolarReflectionManager->iRecPtNum) {
    1145         188 :                 ReflSkySolObs(state.dataSolarReflectionManager->iRecPtNum) = 0.0;
    1146         188 :                 ReflSkySolGnd(state.dataSolarReflectionManager->iRecPtNum) = 0.0;
    1147       17108 :                 for (state.dataSolarReflectionManager->iRayNum = 1;
    1148       34216 :                      state.dataSolarReflectionManager->iRayNum <=
    1149       17108 :                      state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum).NumReflRays;
    1150       16920 :                      ++state.dataSolarReflectionManager->iRayNum) {
    1151       16920 :                     state.dataSolarReflectionManager->HitPntSurfNum =
    1152       16920 :                         state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum)
    1153       16920 :                             .HitPtSurfNum(state.dataSolarReflectionManager->iRayNum, state.dataSolarReflectionManager->iRecPtNum);
    1154             :                     // Skip rays that do not hit an obstruction or ground.
    1155             :                     // (Note that if a downgoing ray does not hit an obstruction it will have HitPtSurfNum = 0
    1156             :                     // if the receiving point is below ground level (see subr. InitSolReflRecSurf); this means
    1157             :                     // that a below-ground-level receiving point receives no ground-reflected radiation although
    1158             :                     // it is allowed to receive obstruction-reflected solar radiation and direct (unreflected)
    1159             :                     // beam and sky solar radiation. As far as reflected solar is concerned, the program does
    1160             :                     // not handle a sloped ground plane or a horizontal ground plane whose level is different
    1161             :                     // from one side of the building to another.)
    1162       16920 :                     if (state.dataSolarReflectionManager->HitPntSurfNum == 0)
    1163        8846 :                         continue; // Ray hits sky or obstruction with receiving pt. below ground level
    1164       16148 :                     state.dataSolarReflectionManager->HitPntRefl =
    1165        8074 :                         state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum)
    1166       16148 :                             .HitPt(state.dataSolarReflectionManager->iRayNum, state.dataSolarReflectionManager->iRecPtNum);
    1167        8074 :                     if (state.dataSolarReflectionManager->HitPntSurfNum > 0) {
    1168             :                         // Ray hits an obstruction
    1169             :                         // Skip hit points on daylighting shelves, from which solar reflection is separately calculated
    1170         824 :                         if (state.dataSurface->SurfDaylightingShelfInd(state.dataSolarReflectionManager->HitPntSurfNum) > 0) continue;
    1171             :                         // Reflected radiance at hit point divided by unobstructed sky diffuse horizontal irradiance
    1172         824 :                         state.dataSolarReflectionManager->HitPtSurfNumX = state.dataSolarReflectionManager->HitPntSurfNum;
    1173             :                         // Each shading surface has a "mirror" duplicate surface facing in the opposite direction.
    1174             :                         // The following gets the correct side of a shading surface in order to get the right value
    1175             :                         // of DifShdgRatioIsoSky (the two sides can have different sky shadowing).
    1176         824 :                         if (state.dataSurface->Surface(state.dataSolarReflectionManager->HitPntSurfNum).IsShadowing) {
    1177        1452 :                             if (dot(state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum)
    1178         484 :                                         .RayVec(state.dataSolarReflectionManager->iRayNum),
    1179         968 :                                     state.dataSurface->Surface(state.dataSolarReflectionManager->HitPntSurfNum).OutNormVec) > 0.0) {
    1180           0 :                                 if (state.dataSolarReflectionManager->HitPntSurfNum + 1 < state.dataSurface->TotSurfaces)
    1181           0 :                                     state.dataSolarReflectionManager->HitPtSurfNumX = state.dataSolarReflectionManager->HitPntSurfNum + 1;
    1182           0 :                                 if (state.dataSurface->SurfDaylightingShelfInd(state.dataSolarReflectionManager->HitPtSurfNumX) > 0) continue;
    1183             :                             }
    1184             :                         }
    1185             : 
    1186         824 :                         if (!state.dataSysVars->DetailedSkyDiffuseAlgorithm || !state.dataSurface->ShadingTransmittanceVaries ||
    1187           0 :                             state.dataHeatBal->SolarDistribution == DataHeatBalance::Shadowing::Minimal) {
    1188         824 :                             state.dataSolarReflectionManager->SkyReflSolRadiance =
    1189        1648 :                                 state.dataSurface->Surface(state.dataSolarReflectionManager->HitPtSurfNumX).ViewFactorSky *
    1190        1648 :                                 state.dataSolarShading->SurfDifShdgRatioIsoSky(state.dataSolarReflectionManager->HitPtSurfNumX) *
    1191         824 :                                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum)
    1192         824 :                                     .HitPtSolRefl(state.dataSolarReflectionManager->iRayNum, state.dataSolarReflectionManager->iRecPtNum);
    1193             :                         } else {
    1194           0 :                             state.dataSolarReflectionManager->SkyReflSolRadiance =
    1195           0 :                                 state.dataSurface->Surface(state.dataSolarReflectionManager->HitPtSurfNumX).ViewFactorSky *
    1196           0 :                                 state.dataSolarShading->SurfDifShdgRatioIsoSkyHRTS(1, 1, state.dataSolarReflectionManager->HitPtSurfNumX) *
    1197           0 :                                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum)
    1198           0 :                                     .HitPtSolRefl(state.dataSolarReflectionManager->iRayNum, state.dataSolarReflectionManager->iRecPtNum);
    1199             :                         }
    1200         824 :                         state.dataSolarReflectionManager->dReflSkySol =
    1201        1648 :                             state.dataSolarReflectionManager->SkyReflSolRadiance *
    1202         824 :                             state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum)
    1203        1648 :                                 .dOmegaRay(state.dataSolarReflectionManager->iRayNum) *
    1204         824 :                             state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum)
    1205        1648 :                                 .CosIncAngRay(state.dataSolarReflectionManager->iRayNum) /
    1206             :                             DataGlobalConstants::Pi;
    1207         824 :                         ReflSkySolObs(state.dataSolarReflectionManager->iRecPtNum) += state.dataSolarReflectionManager->dReflSkySol;
    1208             :                     } else {
    1209             :                         // Ray hits ground;
    1210             :                         // Find radiance at hit point due to reflection of sky diffuse reaching
    1211             :                         // ground directly, i.e., without reflecting from obstructions.
    1212             :                         // Send rays upward from hit point and see which ones are unobstructed and so go to sky.
    1213             :                         // Divide hemisphere centered at ground hit point into elements of altitude Phi and
    1214             :                         // azimuth Theta and create upward-going ray unit vector at each Phi,Theta pair.
    1215             :                         // Phi = 0 at the horizon; Phi = Pi/2 at the zenith.
    1216        7250 :                         Real64 dReflSkyGnd = 0.0; // Factor for ground radiance due to direct sky diffuse reflection
    1217             :                         // Altitude loop
    1218       43500 :                         for (int IPhi = 1; IPhi <= (AltAngStepsForSolReflCalc / 2); ++IPhi) {
    1219             :                             // Third component of ray unit vector in (Theta,Phi) direction
    1220       36250 :                             state.dataSolarReflectionManager->URay(3) = sin_Phi[IPhi];
    1221       36250 :                             Real64 dOmega = cos_Phi[IPhi] * DTheta * DPhi; // Solid angle increment (steradians)
    1222             :                             // Cosine of angle of incidence of ray on ground
    1223       36250 :                             Real64 CosIncAngRayToSky = sin_Phi[IPhi]; // Cosine of incidence angle on ground of ray to sky
    1224             :                             // Azimuth loop
    1225      688750 :                             for (int ITheta = 1; ITheta <= 2 * AzimAngStepsForSolReflCalc; ++ITheta) {
    1226      652500 :                                 state.dataSolarReflectionManager->URay.x = cos_Phi[IPhi] * cos_Theta[ITheta];
    1227      652500 :                                 state.dataSolarReflectionManager->URay.y = cos_Phi[IPhi] * sin_Theta[ITheta];
    1228             :                                 // Does this ray hit an obstruction?
    1229      652500 :                                 hitObs = false;
    1230     3641890 :                                 for (int ObsSurfNum : state.dataSurface->AllShadowPossObstrSurfaceList) {
    1231     3160132 :                                     state.dataSolarReflectionManager->iObsSurfNum = ObsSurfNum;
    1232             :                                     // Horizontal roof surfaces cannot be obstructions for rays from ground
    1233     3160132 :                                     if (state.dataSurface->Surface(state.dataSolarReflectionManager->iObsSurfNum).Tilt < 5.0) continue;
    1234     2420359 :                                     if (!state.dataSurface->Surface(state.dataSolarReflectionManager->iObsSurfNum).IsShadowing) {
    1235     4527158 :                                         if (dot(state.dataSolarReflectionManager->URay,
    1236     4527158 :                                                 state.dataSurface->Surface(state.dataSolarReflectionManager->iObsSurfNum).OutNormVec) >= 0.0)
    1237     1064507 :                                             continue;
    1238             :                                         // Special test for vertical surfaces with URay dot OutNormVec < 0; excludes
    1239             :                                         // case where ground hit point is in back of ObsSurfNum
    1240     2398144 :                                         if (state.dataSurface->Surface(state.dataSolarReflectionManager->iObsSurfNum).Tilt > 89.0 &&
    1241     1199072 :                                             state.dataSurface->Surface(state.dataSolarReflectionManager->iObsSurfNum).Tilt < 91.0) {
    1242     2244292 :                                             state.dataSolarReflectionManager->SurfVert =
    1243     2244292 :                                                 state.dataSurface->Surface(state.dataSolarReflectionManager->iObsSurfNum).Vertex(2);
    1244     2244292 :                                             state.dataSolarReflectionManager->SurfVertToGndPt =
    1245     3366438 :                                                 state.dataSolarReflectionManager->HitPntRefl - state.dataSolarReflectionManager->SurfVert;
    1246     2244292 :                                             if (dot(state.dataSolarReflectionManager->SurfVertToGndPt,
    1247     2244292 :                                                     state.dataSurface->Surface(state.dataSolarReflectionManager->iObsSurfNum).OutNormVec) < 0.0)
    1248      687995 :                                                 continue;
    1249             :                                         }
    1250             :                                     }
    1251     2671428 :                                     PierceSurface(state,
    1252      667857 :                                                   state.dataSolarReflectionManager->iObsSurfNum,
    1253      667857 :                                                   state.dataSolarReflectionManager->HitPntRefl,
    1254      667857 :                                                   state.dataSolarReflectionManager->URay,
    1255      667857 :                                                   state.dataSolarReflectionManager->HitPntObs,
    1256             :                                                   hitObs);
    1257      667857 :                                     if (hitObs) break;
    1258             :                                 }
    1259      652500 :                                 if (hitObs) continue; // Obstruction hit
    1260             :                                 // Sky is hit
    1261      481758 :                                 dReflSkyGnd += CosIncAngRayToSky * dOmega / DataGlobalConstants::Pi;
    1262             :                             } // End of azimuth loop
    1263             :                         }     // End of altitude loop
    1264        7250 :                         ReflSkySolGnd(state.dataSolarReflectionManager->iRecPtNum) +=
    1265        7250 :                             dReflSkyGnd *
    1266        7250 :                             state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum)
    1267       14500 :                                 .dOmegaRay(state.dataSolarReflectionManager->iRayNum) *
    1268        7250 :                             state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum)
    1269       14500 :                                 .CosIncAngRay(state.dataSolarReflectionManager->iRayNum) /
    1270             :                             DataGlobalConstants::Pi;
    1271             :                     } // End of check if ray from receiving point hits obstruction or ground
    1272             :                 }     // End of loop over rays from receiving point
    1273             :             }         // End of loop over receiving points
    1274             : 
    1275             :             // Average over receiving points
    1276          47 :             state.dataSurface->SurfReflFacSkySolObs(state.dataSolarReflectionManager->iSurfNum) = 0.0;
    1277          47 :             state.dataSurface->SurfReflFacSkySolGnd(state.dataSolarReflectionManager->iSurfNum) = 0.0;
    1278          47 :             state.dataSolarReflectionManager->iNumRecPts =
    1279          47 :                 state.dataSolarReflectionManager->SolReflRecSurf(state.dataSolarReflectionManager->iRecSurfNum).NumRecPts;
    1280         235 :             for (state.dataSolarReflectionManager->iRecPtNum = 1;
    1281         235 :                  state.dataSolarReflectionManager->iRecPtNum <= state.dataSolarReflectionManager->iNumRecPts;
    1282         188 :                  ++state.dataSolarReflectionManager->iRecPtNum) {
    1283         188 :                 state.dataSurface->SurfReflFacSkySolObs(state.dataSolarReflectionManager->iSurfNum) +=
    1284         188 :                     ReflSkySolObs(state.dataSolarReflectionManager->iRecPtNum);
    1285         188 :                 state.dataSurface->SurfReflFacSkySolGnd(state.dataSolarReflectionManager->iSurfNum) +=
    1286         188 :                     ReflSkySolGnd(state.dataSolarReflectionManager->iRecPtNum);
    1287             :             }
    1288          47 :             state.dataSurface->SurfReflFacSkySolObs(state.dataSolarReflectionManager->iSurfNum) /= state.dataSolarReflectionManager->iNumRecPts;
    1289          47 :             state.dataSurface->SurfReflFacSkySolGnd(state.dataSolarReflectionManager->iSurfNum) /= state.dataSolarReflectionManager->iNumRecPts;
    1290             :             // Do not allow SurfReflFacBmToDiffSolGnd to exceed the surface's unobstructed ground view factor
    1291          47 :             state.dataSurface->SurfReflFacSkySolGnd(state.dataSolarReflectionManager->iSurfNum) =
    1292          47 :                 min(0.5 * (1.0 - state.dataSurface->Surface(state.dataSolarReflectionManager->iSurfNum).CosTilt),
    1293          47 :                     state.dataSurface->SurfReflFacSkySolGnd(state.dataSolarReflectionManager->iSurfNum));
    1294             :             // Note: the above factors are dimensionless; they are equal to
    1295             :             // (W/m2 reflected solar incident on SurfNum)/(W/m2 unobstructed horizontal sky diffuse irradiance)
    1296             :         } // End of loop over receiving surfaces
    1297           9 :     }
    1298             : 
    1299             : } // namespace SolarReflectionManager
    1300             : 
    1301        2313 : } // namespace EnergyPlus

Generated by: LCOV version 1.13