Nodes of different colours represent the following:
Solid arrows point from a parent (sub)module to the submodule which is descended from it. Dashed arrows point from a module being used to the module or program unit using it. Where possible, edges connecting nodes are given different colours to make them easier to distinguish in large graphs.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
type(BSDFRefPoints), | intent(inout) | :: | StateRefPoint | |||
type(BSDFRefPointsGeomDescr), | intent(inout) | :: | DaylghtGeomDescr | |||
integer, | intent(in) | :: | ZoneNum | |||
integer, | intent(in) | :: | iWin | |||
real(kind=r64), | intent(in), | dimension(3) | :: | RefPoint | ||
integer, | intent(in) | :: | curFenState | |||
integer, | intent(in) | :: | NBasis | |||
integer, | intent(in) | :: | NTrnBasis | |||
real(kind=r64), | intent(in) | :: | AZVIEW | |||
integer, | intent(in) | :: | NWX | |||
integer, | intent(in) | :: | NWY | |||
real(kind=r64), | intent(in), | dimension(3) | :: | W2 | ||
real(kind=r64), | intent(in), | dimension(3) | :: | W21 | ||
real(kind=r64), | intent(in), | dimension(3) | :: | W23 | ||
real(kind=r64), | intent(in) | :: | DWX | |||
real(kind=r64), | intent(in) | :: | DWY | |||
real(kind=r64), | intent(in), | dimension(3) | :: | WNorm | ||
real(kind=r64), | intent(in) | :: | WinElArea | |||
integer, | intent(in) | :: | CalledFrom | |||
integer, | intent(in), | optional | :: | MapNum |
Nodes of different colours represent the following:
Solid arrows point from a procedure to one which it calls. Dashed arrows point from an interface to procedures which implement that interface. This could include the module procedures in a generic interface or the implementation in a submodule of an interface in a parent module. Where possible, edges connecting nodes are given different colours to make them easier to distinguish in large graphs.
Nodes of different colours represent the following:
Solid arrows point from a procedure to one which it calls. Dashed arrows point from an interface to procedures which implement that interface. This could include the module procedures in a generic interface or the implementation in a submodule of an interface in a parent module. Where possible, edges connecting nodes are given different colours to make them easier to distinguish in large graphs.
SUBROUTINE InitializeCFSStateData(StateRefPoint, DaylghtGeomDescr, ZoneNum, iWin, RefPoint, curFenState, NBasis, NTrnBasis, &
AZVIEW, NWX, NWY, W2, W21, W23, DWX, DWY, WNorm, WinElArea, CalledFrom, MapNum)
! SUBROUTINE INFORMATION:
! AUTHOR Simon Vidanovic
! DATE WRITTEN June 2013
! MODIFIED na
! RE-ENGINEERED na
! PURPOSE OF THIS SUBROUTINE:
! Initialize daylight state data for current
! METHODOLOGY EMPLOYED:
! <description>
! REFERENCES:
! na
! USE STATEMENTS:
use vectors
use WindowComplexManager, only: PierceSurfaceVector, DaylghtAltAndAzimuth
implicit none ! Enforce explicit typing of all variables in this routine
! SUBROUTINE ARGUMENT DEFINITIONS:
type(BSDFRefPoints), intent(inout) :: StateRefPoint
type(BSDFRefPointsGeomDescr), intent(inout) :: DaylghtGeomDescr
integer, intent(in) :: ZoneNum ! Current zone number
integer, intent(in) :: iWin
real(r64), dimension(3), intent(in) :: RefPoint ! reference point
integer, intent(in) :: curFenState
integer, intent(in) :: NBasis
integer, intent(in) :: NTrnBasis
!integer, intent(in) :: NRefPt
integer, intent(in) :: CalledFrom
real(r64), intent(in) :: AZVIEW
integer, intent(in) :: NWX
integer, intent(in) :: NWY
real(r64), intent(in), dimension(3) :: W2
real(r64), intent(in), dimension(3) :: W21
real(r64), intent(in), dimension(3) :: W23
real(r64), intent(in) :: DWX
real(r64), intent(in) :: DWY
real(r64), intent(in), dimension(3) :: WNorm ! unit vector from window (point towards outside)
real(r64), intent(in) :: WinElArea
integer, intent(in), optional :: MapNum
! SUBROUTINE LOCAL VARIABLES
type(vector) :: Centroid ! current window element centroid
integer :: curWinEl
integer :: IRay
integer :: iHit
integer :: TotHits
integer :: JSurf
real(r64) :: DotProd !Temporary variable for manipulating dot product .dot.
integer :: NSky
integer :: NGnd
integer :: NReflSurf
integer :: MaxTotHits
integer :: IX
integer :: IY
real(r64), dimension(3) :: RWin ! window element center point (same as centroid)
type(vector) :: HitPt ! surface hit point
real(r64), dimension(3) :: V !vector array
real(r64) :: LeastHitDsq ! dist^2 from window element center to hit point
real(r64) :: HitDsq
real(r64) :: TransRSurf
integer :: I, J
real(r64), dimension(3) :: GroundHitPt ! Coordinates of point that ray hits ground (m)
! Refrence point data
! integer :: iRefPoint
! temporary arrays for surfaces
integer, dimension(:), allocatable :: TmpSkyInd !Temporary sky index list
integer, dimension(:), allocatable :: TmpGndInd !Temporary gnd index list
type (vector), dimension(:), allocatable :: TmpGndPt !Temporary ground intersection list
real(r64), dimension(:), allocatable :: TmpGndMultiplier !Temporary ground obstruction multiplier
integer, dimension(:), allocatable :: TmpRfSfInd !Temporary RefSurfIndex
integer, dimension(:), allocatable :: TmpRfRyNH !Temporary RefRayNHits
integer, dimension(:,:), allocatable :: TmpHSurfNo !Temporary HitSurfNo
type(vector), dimension(:,:), allocatable :: TmpHitPt !Temporary HitPt
real(r64), dimension(:,:), allocatable :: TmpHSurfDSq !Temporary HitSurfDSq
! each complex fenestration state can have different number of basis elements. This is the reason for making
! temporary arrays allocated at this moment. It is also important to note that deallocation needs to occur at the
! end of current fenestration state loop
if (.not.allocated(TmpSkyInd)) allocate(TmpSkyInd(NBasis))
TmpSkyInd = 0.0d0
if (.not.allocated(TmpGndInd)) allocate(TmpGndInd(NBasis))
TmpGndInd = 0.0d0
if (.not.allocated(TmpGndPt)) allocate(TmpGndPt(NBasis))
TmpGndPt = vector(0.0d0, 0.0d0, 0.0d0)
if (.not.allocated(TmpGndMultiplier)) allocate(TmpGndMultiplier(NBasis))
TmpGndMultiplier = 0.0d0
if (.not.allocated(TmpRfSfInd)) allocate(TmpRfSfInd(NBasis))
TmpRfSfInd = 0.0d0
if (.not.allocated(TmpRfRyNH)) allocate(TmpRfRyNH(NBasis))
TmpRfRyNH = 0.0d0
if (.not.allocated(TmpHSurfNo)) allocate(TmpHSurfNo(NBasis,TotSurfaces))
TmpHSurfNo = 0.0d0
if (.not.allocated(TmpHitPt)) allocate(TmpHitPt(NBasis,TotSurfaces))
TmpHitPt = vector(0.0d0, 0.0d0, 0.0d0)
if (.not.allocated(TmpHSurfDSq)) allocate(TmpHSurfDSq(NBasis,TotSurfaces))
TmpHSurfDSq = 0.0d0
! find if reference point belongs to light tube of outgoing bsdf direction. This works for entire window and not window
! elements.
! initialization for each reference point
!do iRefPoint = 1, NRefPt
!if (CalledFrom == CalledForRefPoint) then
! RefPoint = ZoneDaylight(ZoneNum)%DaylRefPtAbsCoord(iRefPoint, 1:3)
!else
! RefPoint = IllumMapCalc(MapNum)%MapRefPtAbsCoord(irefPoint, 1:3)
!end if
call CFSRefPointPosFactor(RefPoint, StateRefPoint, iWin, CurFenState, NTrnBasis, AZVIEW)
!end do
curWinEl = 0
! loop through window elements. This will calculate sky, ground and reflection bins for each window element
do IX = 1, NWX
do IY = 1, NWY
curWinEl = curWinEl + 1
! centroid coordinates for current window element
Centroid = W2 + (REAL(IX,r64) - 0.5d0) * W23 * DWX + (REAL(IY,r64) - 0.5d0) * W21 * DWY
RWin = Centroid
!do iRefPoint = 1, NRefPt
!if (CalledFrom == CalledForRefPoint) then
! RefPoint = ZoneDaylight(ZoneNum)%DaylRefPtAbsCoord(iRefPoint, 1:3)
!else
! RefPoint = IllumMapCalc(MapNum)%MapRefPtAbsCoord(iRefPoint, 1:3)
!end if
call CFSRefPointSolidAngle(RefPoint, RWin, WNorm, StateRefPoint, &
DaylghtGeomDescr, iWin, CurFenState, &
NTrnBasis, curWinEl, WinElArea)
!end do
NSky = 0
NGnd = 0
NReflSurf = 0
MaxTotHits = 0
! Calculation of potential surface obstruction for each incoming direction
do IRay = 1, NBasis
iHit = 0
TotHits = 0
do JSurf = 1, TotSurfaces
! the following test will cycle on anything except exterior surfaces and shading surfaces
if( Surface(JSurf)%HeatTransSurf .AND. Surface(JSurf)%ExtBoundCond /= ExternalEnvironment) cycle
! skip the base surface containing the window and any other subsurfaces of that surface
if( JSurf == Surface(IWin)%BaseSurf .OR. Surface(JSurf)%BaseSurf == Surface(IWin)%BaseSurf) cycle
! skip surfaces that face away from the window
DotProd = ComplexWind(IWin)%Geom(CurFenState)%sInc(IRay) .dot. Surface(JSurf)%NewellSurfaceNormalVector
if( DotProd >= 0 ) cycle
call PierceSurfaceVector(JSurf, Centroid, ComplexWind(IWin)%Geom(CurFenState)%sInc(IRay), IHit, HitPt)
if (IHit <= 0) cycle
IHit = 0 !A hit, clear the hit flag for the next cycle
if (TotHits == 0) then
! First hit for this ray
TotHits = 1
NReflSurf = NReflSurf + 1
TmpRfSfInd(NReflSurf) = IRay
TmpRfRyNH(NReflSurf) = 1
TmpHSurfNo(NReflSurf,1) = JSurf
TmpHitPt(NReflSurf,1) = HitPt
V = HitPt - Centroid !vector array from window ctr to hit pt
LeastHitDsq = DOT_PRODUCT (V , V) !dist^2 window ctr to hit pt
TmpHSurfDSq(NReflSurf , 1) = LeastHitDsq
if(.NOT.Surface(JSurf)%HeatTransSurf .AND. Surface(JSurf)%SchedShadowSurfIndex /= 0) then
TransRSurf = 1.0d0 !If a shadowing surface may have a scheduled transmittance,
! treat it here as completely transparent
else
TransRSurf= 0.0d0
end if
else
V = HitPt - Centroid
HitDsq = DOT_PRODUCT (V , V)
if (HitDsq >= LeastHitDsq) then
if (TransRSurf > 0.0d0) then !forget the new hit if the closer hit is opaque
J = TotHits + 1
if (TotHits > 1) then
do I = 2, TotHits
if ( HitDsq < TmpHSurfDSq(NReflSurf , I) ) THEN
J = I
exit
end if
end do
if (.NOT.Surface(JSurf)%HeatTransSurf .AND. &
Surface(JSurf)%SchedShadowSurfIndex == 0) then
! The new hit is opaque, so we can drop all the hits further away
TmpHSurfNo (NReflSurf , J) = JSurf
TmpHitPt (NReflSurf , J) = HitPt
TmpHSurfDSq (NReflSurf , J) = HitDsq
TotHits = J
else
! The new hit is scheduled (presumed transparent), so keep the more distant hits
! Note that all the hists in the list will be transparent except the last,
! which may be either transparent or opaque
if (TotHits >= J ) then
do I = TotHits , J , -1
TmpHSurfNo (NReflSurf , I+1) = TmpHSurfNo (NReflSurf , I)
TmpHitPt (NReflSurf , I+1) = TmpHitPt (NReflSurf , I)
TmpHSurfDSq (NReflSurf , I+1) = TmpHSurfDSq (NReflSurf , I)
end do
TmpHSurfNo (NReflSurf , J) = JSurf
TmpHitPt (NReflSurf , J) = HitPt
TmpHSurfDSq (NReflSurf , J) = HitDsq
TotHits = TotHits + 1
end if
end if ! if (.NOT.Surface(JSurf)%HeatTransSurf .AND. Surface(JSurf)%SchedShadowSurfIndex == 0) then
end if ! if (TotHits > 1) then
end if ! if (TransRSurf > 0.0d0) then
else ! if (HitDsq >= LeastHitDsq) then
! A new closest hit. If it is opaque, drop the current hit list,
! otherwise add it at the front
LeastHitDsq = HitDsq
if(.NOT.Surface(JSurf)%HeatTransSurf .AND. Surface(JSurf)&
%SchedShadowSurfIndex /= 0) then
TransRSurf = 1.0d0 ! New closest hit is transparent, keep the existing hit list
do I = TotHits , 1 , -1
TmpHSurfNo (NReflSurf , I+1) = TmpHSurfNo (NReflSurf , I)
TmpHitPt (NReflSurf , I+1) = TmpHitPt (NReflSurf , I)
TmpHSurfDSq (NReflSurf , I+1) = TmpHSurfDSq (NReflSurf , I)
TotHits = TotHits + 1
end do
else
TransRSurf = 0.0d0 ! New closest hit is opaque, drop the existing hit list
TotHits = 1
end if
TmpHSurfNo (NReflSurf,1) = JSurf ! In either case the new hit is put in position 1
TmpHitPt (NReflSurf,1) = HitPt
TmpHSurfDSq (NReflSurf, 1) = LeastHitDsq
end if
end if
end do ! do JSurf = 1, TotSurfaces
if (TotHits <= 0 ) then
!This ray reached the sky or ground unobstructed
if (ComplexWind(IWin)%Geom(CurFenState)%sInc(IRay)%z < 0.0d0) then
!A ground ray
NGnd = NGnd + 1
TmpGndInd(NGnd) = IRay
TmpGndPt(NGnd)%x = Centroid%x - (ComplexWind(IWin)%Geom(CurFenState)%sInc(IRay)%x / &
ComplexWind(IWin)%Geom(CurFenState)%sInc(IRay)%z) * Centroid%z
TmpGndPt(NGnd)%y = Centroid%y - (ComplexWind(IWin)%Geom(CurFenState)%sInc(IRay)%y / &
ComplexWind(IWin)%Geom(CurFenState)%sInc(IRay)%z) * Centroid%z
TmpGndPt(NGnd)%z = 0.0d0
! for solar reflectance calculations, need to precalculate obstruction multipliers
if (CalcSolRefl) then
GroundHitPt = TmpGndPt(NGnd)
TmpGndMultiplier(NGnd) = CalcObstrMultiplier(GroundHitPt, AltAngStepsForSolReflCalc, AzimAngStepsForSolReflCalc)
end if
else
!A sky ray
NSky = NSky + 1
TmpSkyInd(NSky) = IRay
end if
else
!Save the number of hits for this ray
TmpRfRyNH (NReflSurf) = TotHits
end if
MaxTotHits = MAX(MaxTotHits, TotHits)
end do ! do IRay = 1, ComplexWind(IWin)%Geom(CurFenState)%Inc%NBasis
! Fill up state data for current window element data
StateRefPoint%NSky(curWinEl) = NSky
StateRefPoint%SkyIndex(curWinEl, 1:NSky) = TmpSkyInd(1:NSky)
StateRefPoint%NGnd(curWinEl) = NGnd
StateRefPoint%GndIndex(curWinEl, 1:NGnd) = TmpGndInd(1:NGnd)
StateRefPoint%GndPt(curWinEl, 1:NGnd) = TmpGndPt(1:NGnd)
StateRefPoint%GndObstrMultiplier(curWinEl, 1:NGnd) = TmpGndMultiplier(1:NGnd)
StateRefPoint%NReflSurf(curWinEl) = NReflSurf
StateRefPoint%RefSurfIndex(curWinEl, 1:NReflSurf) = TmpRfSfInd(1:NReflSurf)
StateRefPoint%RefRayNHits(curWinEl, 1:NReflSurf) = TmpRfRyNH(1:NReflSurf)
StateRefPoint%HitSurfNo(curWinEl, 1:NReflSurf, 1:MaxTotHits) = TmpHSurfNo(1:NReflSurf, 1:MaxTotHits)
StateRefPoint%HitSurfDSq(curWinEl, 1:NReflSurf, 1:MaxTotHits) = TmpHSurfDSq(1:NReflSurf, 1:MaxTotHits)
StateRefPoint%HitPt(curWinEl, 1:NReflSurf, 1:MaxTotHits) = TmpHitPt(1:NReflSurf, 1:MaxTotHits)
end do ! do IY = 1, NWY
end do ! do IX = 1, NWX
! Deallocate temporary arrays
if (allocated(TmpSkyInd)) deallocate(TmpSkyInd)
if (allocated(TmpGndInd)) deallocate(TmpGndInd)
if (allocated(TmpGndPt)) deallocate(TmpGndPt)
if (allocated(TmpGndMultiplier)) deallocate(TmpGndMultiplier)
if (allocated(TmpRfSfInd)) deallocate(TmpRfSfInd)
if (allocated(TmpRfRyNH)) deallocate(TmpRfRyNH)
if (allocated(TmpHSurfNo)) deallocate(TmpHSurfNo)
if (allocated(TmpHitPt)) deallocate(TmpHitPt)
if (allocated(TmpHSurfDSq)) deallocate(TmpHSurfDSq)
END SUBROUTINE InitializeCFSStateData