# Optical Losses

These routines compute all optical losses due to shading, soiling, snow, and refraction through the module’s glass superstrate. These losses should be computed in the order that they physically occur so that the energy differences and annualized loss ratios can be correctly accounted, namely:

• Lost energy due far (horizon or ridgeline) shading, applied to the direct beam component
• Lost energy due to near row-on-row direct beam shading (obstructed direct beam light that never reaches the module)
• Lost energy due to diffuse shading (obstructed diffuse light due to nearby objects)
• Lost energy due to snow and soiling (energy incident on the panel that never reaches the glass surface
• Lost energy due to IAM (energy incident on the glass that never reaches the active semiconductor surface)
• Lost energy due to spectral mismatch between the module response to the G173 spectrum and the actual spectrum

## Algorithm

Compute the effective irradiance as seen by the semiconductor, after all of the shading, reflection, refraction, and soiling losses have been subtracted from the direct beam, diffuse sky, and diffuse albedo (ground reflected) irradiance components.

1.) If the apparent solar elevation angle is lower than the linearly-interpolated horizon line, for the current solar azimuth, set the plane-of-array beam irradiance component to zero. This is a binary operation affecting only the beam component. The diffuse components are untouched.

(This is equivalent to saying that if the zenith angle is greater than the compliment of the horizon elevation angle, then the effective beam irradiance on the plane is zero.)

2.) Compute the effective direct beam irradiance on the tilted plane after the incidence angle modifier (IAM) and row-on-row shading losses.

3.) Compute the effective diffuse sky irradiance on the tilted plane after diffuse IAM and shading losses.

4.) Compute the effective diffuse ground irradiance on the tilted plane after diffuse IAM and shading

5.) Find the total effective irradiance after the effect of soiling, snow, and other contaminants.

6.) Find the total effective irradiance available for PV conversion, after spectral mismatch.

These series of function compute the shading factors on the different irradiance components, for unlimited PV sheds. PlantPredict supports various shading and incident angle modifier algorithms.

“3D” is the short hand way of referring to the depth pass stencil buffer shadow volume (DPSBSV) algorithm used for advanced shadow modeling in PlantPredict. It is carried out with the aid of the OpenGL API. The DPSBSV is one of the most common ways that shadows are calculated in video games today. It is a graphical solution, meaning the algorithm must first draw out all of the objects in a given DC field (i.e. all PV tables, as well as any external objects added which could cause shading), and after it has done so it will count up the number of pixels that are in shadow. This is carried out in a 5 step process.

1.) Create the “shadow volumes”

Shadow volumes are objects whose interior represents where the shadow is cast. This is done by projecting the vertices of a solar panel or object in the opposite direction to the sun position vector. Figure X: Objects without shadows

Figure Xa&b: Shadows rendered as solid shadow volume

2.) Render the solar panels

3.) Render the shadow volume into the stencil buffer

The stencil buffer allows the user to make decisions on which areas to draw. Every stencil buffer pixel starts at 0; if a pixel sees the top of a shadow, the stencil buffer is incremented. If a pixel sees the bottom of a shadow, the stencil buffer is decremented. If it only sees a solar panel it is kept at 0.

4.) Redraw panels everywhere the stencil buffer is equal to 0

5.) Run occlusion test to get linear shading factor

So that it can perform the shading calculations as quickly as possible, the algorithm performs what is known as an occlusion test. This purpose of the test is to determine which hours have shading and to exclude any hours without shading from the calculation. Once an hour passes the occlusion test the linear shading factor can be determined simply by taking a ratio of the sum of the stencil buffer (shown in step 4) to the original field rendered with no shadows (shown in step 2).

3D plus Electrical Shading

For electrical shading the operation is a bit more complex. Instead of just taking a simple ratio of shaded area, the algorithm must determine which electrical strings are shaded. The number of total electrical strings in a DC Field is known based on the inputs, but In order to uniquely identify all electrical strings that are shaded the algorithm renders each one as a different color and passes these colors (represented numerically) into the stencil buffer. Finally, the algorithm counts up the total number of unique colors and divides it by the total number of electrical strings to yield the electrical shading factor for a given hour.

At this current point in time, the 3D shading algorithm is only calculating the shading factors for the direct beam component of irradiance. The impact on diffuse irradiance is still being calculated by PlantPredict’s original diffuse model which accounts for the impact of row-to-row shading and slope. The PlantPredict team is working to quantify an accurate impact of diffuse shading in the 3D model, and how it pertains to the addition of near field objects.

Application of Shading Factors

Linear –

The hourly linear shading factors are multiplied by the incident beam irradiance values so that a linear effective beam irradiance value can be determined.

Electrical –

The hourly electrical shading factors are multiplied by the incident beam irradiance values so that an electrical effective beam irradiance value can be determined.

# Direct Beam Horizon Shading

## Algorithm Introduction

Satellite-derived irradiance data usually does not account for local topography. Horizon shading occurs when the beam component of the irradiance is occluded by a ridgeline, usually near the times of sunrise or sunset. The apparent sunrise and sunset times may be delayed and advanced, respectively, if mountains are present.

The shading effect is treated globally for the entire power plant, in other words, partial shading of a plant by distant mountains is not taken into account. This algorithm outputs a scaling value (between 0 and 1) that factors the direct beam irradiance proportional to the amount of time the sun was occluded during that simulation interval.

To estimate the amount of time the sun is blocked from far shading we use the solar geometry, the horizon profile (hProfile) and linear interpolation at 1-minute resolutions.

### The following summarizes the algorithm steps:

1.) Determine if the current timestamp could be subject to horizon shading

Using linear interpolation we can find the ridgeline elevation at any azimuth. Comparing these elevations with the sun’s elevations (at the same azimuth) over the course of the entire interval, we can determine if any horizon shading occurs.

The first step in this algorithm is to compile a list of the solar azimuth-elevation pairs (SEtSAt). If lowest solar elevation angle (has to be greater than 0) of the previous, current and next timestep is smaller than the highest elevation in the ridgeline profile then we will run the horizon shading analysis, otherwise we assume there is no shading.

2.) Estimate the time that the sun traverses the ridgeline

To determine if and when the sun traverses the ridgeline, we need to start with the timestamp that has the lowest solar elevation. The timestamps with the lowest solar elevation will always occur at the start (sunrise) or end (sunset) of a given interval. From the rotation of the earth, we know:

If 0 ≤ SAt ≤ 180 then the sun is rising and the lowest elevation is at the start of the interval, otherwise it’s at the end of the interval. For example, in the figure below, the sun is rising and the current timestamp is 7:30 AM, therefore the lowest elevation of the interval occurs at 7:00 AM. Sometime the lowest elevation of the interval is the Sunrise or Sunset time, this occurs when a sunrise or sunset event occurs during the current simulation interval.

With starting time, we can interpolate to find the solar azimuth and elevation at that time, and in turn the azimuth can be used to find the ridgeline elevation. While the ridgeline elevation is greater the solar elevation, the sun is still being obscured and we can step (interpolate) the next minute’s solar azimuth-elevation pair. Once we find the minute where the sun traverses the elevation profile, we store that time and break from the loop.

### Figure 9. Example for Interpolating the Apparent Sunrise Time

3.) Estimate amount of beam irradiance that was obscured by the ridgeline

Now that we know the solar transit time above the ridgeline time, we need to calculate how much beam irradiance was blocked by the ridgeline. Consider the configuration diagramed in figure 9 above. Our current middle-of-irradiance interval timestamp is 7:30 AM and the transit time is 7:15 AM. If we found the transit time to be before 7:00 AM then we would have no shading to calculate, but in this case we know that the sun was obscured for 15 minutes of the current interval. We then estimate the amount of beam shading by finding what fraction of the interval was occluded, in this case 15 min / 60 min or 25%. This corresponds to a beam shading factor or SFactor of 0.72.

There are shortcomings of this algorithm, namely that it tends to overpredict beam shading because it treats the irradiance over every minute of interval as equal. In reality, the irradiance varies throughout the hour and therefor the shading factor should reflect this. This could be achieved in the future by interpolating the irradiances at t-1,t and t+1 to come up with a more accurate shading factor value.

Additionally, the algorithm could be improved by reducing the diffuse ground and diffuse sky irradiances but both of these changes need to be motivated by measured data before they can be implemented.

# Electrical Shading Models and Calculation

## Fractional Effect –

This model assumes that all module strings in a table are strung together according to vertical rank. Below, an example portrays a PV table that consists of 2 strings each consisting of 8 modules.

The table is partially shaded on the bottom corner, and due to the fact that at least one module in the bottom string gets shaded the entire string experiences power loss due to electrical shading. The fractional effect value determines how much electrical shading power loss will be applied. A fractional effect of 100% will have the entire strings power go to zero once a portion (of greater than 1%) gets shaded. This will determine the electrical beam shading factor that gets calculated for the DC field and, in turn, the effective irradiance that the modules will see (explained further in Shade Factor Calculation).

A 75% fractional effect will reduce the entire string’s electrical beam shading factor to 25% of the linear beam shading factor, essentially reducing the loss due to electrical shading by 25%.

# Shade Factor Calculation:

The below figure defines certain variables used in the Shade Factor (UshBES) calculation. Xm and Ym are the module’s length and width dimension, respectively. The other two dimensions, Xs and Ys represent the distance, length and width respectively, that the shadow overlaps the module.

Knowing these dimensions at each time step, PlantPredict is able to calculate the shading factor using the following equations:

In which Z is the shading response table (examples of two such tables were provided in the c-SI 3 Diode Default and Custom Shading Table sections above). Once the shading factor has been determined, PlantPredict will then use that to compute the Global Effective Beam Irradiance (Ee,ES) after it has been adjusted for electrical shading.

UshBES,filt is a modified version of the shading factor values (UshBES) that were calculated earlier, but filtered to convert any value that is to 1 (no shading applied) if the area of the module shaded is less than 1%.

The next tasks that the logic performs is to recalculate the cell temperature to reflect the fact that the module is being subjected to a new Global Effective Beam Irradiance value. This calculation is done based on the Module Temperature model that is selected in the simulation settings. To see this calculation please refer to the Module Temperature section lower in this document.

Once Ee,ES and Tcell,mod are determined, the last step is for PlantPredict to plug these into the functions that compute the 5 parameters of the I-V curve and subsequently the single diode model so that the maximum power point energy after electrical shading (EMPP,ES) can be found.

The Electrical Loss % (as seen in the results screen) is calculated as follows:

In which EMPP is the maximum power point energy if it were to be calculated with linear shading.

# Near Direct Beam Shading

This series of functions computes the shading factors on the direct beam component of irradiance. Given the sun’s azimuth angle and zenith angle (elevation angle’s complement), as well as the array’s orientation, spacing, tilt, and dimensions, the functions compute the geometry of the shadow cast by one row on an adjacent row as well as the fraction of the array that is shaded by tables within the array. Edge effects can be excluded based on user selection, which represents ignoring the unshaded edges of rows for DC fields that do not represent full arrays, typically used for the case when the array is designed for string inverters.

Figure 10 illustrates the direct beam row-to-row shading. The current tilt of the modules is represented here by 𝛽 for the fixed tilt or tracker tilt angle.

### Figure 10. One-dimensional shed shading illustration of either a fixed-tilt or tracker, as viewed from the edge of the row

The following diagram, Figure 11, illustrates the variables used to calculate the direct beam shading on an array. In the figure, Δpp is the distance between tables, LT is the length of each table, and WT is the width of each table. Nr is the number of sheds. Although not shown, module azimuth angle, γ, is defined with due north being 0°, and is dependent upon user input DC field orientation. Also not shown, NT is the number of tables per shed. Although Figure 11 illustrates a fixed tilt array, all geometric definitions are the same for horizontal tracking arrays.

### Figure 11. Two-dimensional mutual shed-on-shed shading. Example shows a fixed-tilt array with 180 azimuth (south-facing). Diagram adapted from: Dan Weinstock, Joseph Appelbaum, Shadow Variation on Photovoltaic Collectors in a Solar Field. Proceedings of the 23rd IEEE Convention of Electrical and Electronics Engineers in Israel, 2004.

Figure 12 shows the geometries that define a terrain slope. The angle of the slope is defined as 𝛼. This can be measured as the angle between the vector normal to the surface of the terrain and the zenith. The direction faced by the slope, called the aspect, is the azimuth from north of the vector normal to the surface of the terrain, γg.

### Figure 12. Geometry of a slope. Terrain slope 𝛼 and aspect γg. Figure adapted and recreated from http://www.ces.iisc.ernet.in/grass/grass72/manuals/i.topo.corr.html

The diagram in Figure 13 illustrates the additional variables of slope 𝛼 and aspect γg used for modeling an array that is built along a slope. The module azimuth of the array points south, but the slope geometry is independent of the module azimuth.

## Inputs

[1] The incidence angle is not used to compute the shading factor per se, but is a test to verify if the sun is behind the plane of the array.

## Algorithm (Vector Method)

### Decomposition of a Slope

The terrain slope angle 𝛼 can be defined by the user for any aspect, or azimuth direction, γg. The aspect is defined according to the same convention as module azimuth angle.

To evaluate the shading that results from such a slope, it is helpful to decompose the slope into the two separate slopes that are along and perpendicular to the azimuth angle of the modules.

The component of the slope that is along the module azimuth angle affects the shading that occurs from one row onto another, or row-to-row shading. Since this shading occurs across rows, it is being termed Transverse shading. The contours of the slope in the transverse direction run parallel to the rows.

The component of the slope perpendicular to the module azimuth angle is along the rows of the array. It causes the tables of the rows to be at staggered heights, for a horizontal tracker system. This slope therefore creates shading between neighboring tables along a row, which is being termed Adjacent shading, or table-to-table shading. The overall shading impact of the slope 𝛼 is the shading due to both components. Depending on the relative location of a given table in the array, it might cast either Transverse shading, Adjacent shading, or a combination of the two.

The component of the slope that affects the Transverse shading between neighboring rows at a positive tilt is defined as transverse slope, 𝜓. Transverse shading is illustrated in Figure 14.

The component of the slope that affects the Adjacent shading between neighboring tables is defined as adjacent slope, 𝜙. For a tracker system, an adjacent slope facing away from the equator is considered positive. For a fixed tilt system, a westward facing adjacent slope is considered positive. Figure 15 illustrates the adjacent shading created between tables within a row at different heights along the slope.

### Figure 15. Treatment of sheds on terrain with an adjacent slope, 𝜙, showing the horizontal table-to-table spacing Lg of the slope angle 𝜓. The figure could be interpreted as either: (a) the front-view (looking east) of a northern-hemisphere horizontal tracker row on a north-facing slope at noon (sun due south at right), or alternatively: (b) looking north at a fixed tilt row on an west-facing slope at sunrise.

The components of the slope are decomposed into the perpendicular directions with the following calculations. The ground slope α can be positive or negative, while the aspect γg must be positive, from 0 to 2𝜋.

For a fixed tilt system with a positive latitude:

For a fixed tilt system with a negative latitude or a tracker system with a positive latitude:

For a tracker system with a negative latitude:

The generalized shading fraction is calculated in the following algorithm.

1. First, define a right-handed table-centric Cartesian coordinate system for a 2-row geometry:

Origin:    Center of second table, at axis of rotation.

i:              positive in the direction of the first table (i.e. positive in the south direction for south-facing arrays).

j:              positive to the right along the lengthwise centerline of second table (from perspective of observer facing table).

k:             positive in the direction of the zenith.

The table causing shading in the row-to-row direction or the table-to-table direction (lengthwise along the rows) is the primary table, PTp or Prp, respectively, or Pp in general. The table that is shaded is the secondary table, Sp.

Figure 16. Primary tables Pp relative to the shaded table Sp at a time other than solar noon. Primary table across the rows (or within the column), Prp, which creates Transverse shading. Primary table along (within) the row, PTp , which causes Adjacent (table-to-table) shading. The figure could be interpreted as either: (a) a N-S axis tracker system, or alternatively: (b) a south-facing fixed tilt system.

2. The solar unit vector s, transformed from spherical to the table-centric Cartesian coordinates, is:

3. To calculate the shadow, we need to determine the vertical spacing between a corresponding point on the primary table (Pp) and the shaded table (Sp), conveniently determined by the height difference at the tops of corresponding posts of the tables. This spacing is defined as the glass-to-glass height difference Δgg.

4. The general equation of a plane is defined as a*x + b*y + c*z = 0. In our case, the tilted PV plane Sp intersects the plane normal to the surface of the earth at the y axis of our coordinate system, so the second row’s table sits in the plane a*x + c*z = 0, where a = sin(𝛽) and c = cos(𝛽). The equation for the plane in which the secondary table lies is therefore:

### Figure 17. Detailed shaded table  geometry and coordinate system.

5. To find the percentage shaded as a function of the solar unit vector, first define the four vertices (the four corners) of the primary (Pp) table.

### Figure 18. Vertex assignments of primary (Pp) table

The i,j,k coordinates of these vertices are given in the following matrix. The vertices are listed in the rows of the matrix, starting on the top-left of the table, and proceeding clockwise.

6. Find the projection of the solar unit vector at each of these points Pp onto the plane containing the secondary table amounts to findthe intersection of the plane containing the secondary table with solar rays grazing each vertex of the primary table.

### Figure 19. The projection of each of the points Prp along the solar unit vector s onto the plane containing the secondary table, where r is distance the solar ray has to travel between the two planes. This is at a time other than solar noon.

This can be evaluated as follows. The parameterization of the line made up of the solar ray passing through each vertex, n, as a function of r is given by:

Here s is the solar unit vector. Substitute the formula for the line into the equation of the plane:

Here r has become the distance the solar ray has to travel between the two tables. Determine the value of r, where a = sin 𝛽 and b = cos 𝛽:

Substitute r into the formula for the line to find the point of intersection in the secondary plane of the solar rays through each of the primary table’s vertices.

7. Project the shading points Sp from 3D to 2D space, where y runs along the length of the table and x runs along the width of the table, as shown in Figure 20.

### Figure 20. Relative locations of the shadow vertices P (shadow greyed) on the shaded table in 2D space.

In this coordinate system, y is simply the y-coordinate in the 3D panel-centric coordinate system, plus half the table length; x is the hypotenuse of i and k, added or subtracted from half the table width depending on the direction of i relative to the center of the table:

8. If only one shadow falls on the table, it would be sufficient to determine the location of the corner of the shadow within the table to calculate the portion of the table that is shaded. However, if there might be multiple shadows on the table, the location of each shadow must be determined.

Next, find the location of the corners of PTn, the rectangular subset of the shadow that falls within the secondary table. The x-component of the shadow that exists on the secondary table is determined by testing the x-component of the first and third vertices of the shadow P against the x-component of the first and third vertices of the secondary table. The distance y is similarly determined with the same vertices. PTn can be calculated using either pair of opposite corners. It is assumed the planes of the tables are parallel to each other. The 2-dimensional coordinates of the secondary table’s vertices are given by V.

### Figure 21. Two-dimensional row shading on the table (shadow greyed)

9. The linear shading fraction SF of a single shadow on one table is then found by2:

[2] Finding the shaded fraction could be implemented much more simply by using MATLAB’s rectangle intersection function Rectint, which requires as its arguments only the vertex vectors of polygons Sp and Pp.

10. For an entire array comprised of Nr rows or Nt columns (N in general) causing mutual shading, spaced equally, the shading fraction is slightly less because the first row or column remains unshaded, while the others are affected equally by their mutual shading.

If multiple shadows from additional tables are considered, then μsh,B is calculated by a weighted average of the shading fraction on the relative number of tables.

For modeling portions of a larger array that has varying DC Field configurations, or when multiple small inverters and corresponding DC Fields make up an array, it may be preferred to assume that no unshaded columns- tables along the edge of the DC Field- exist. This is accomplished by assuming the Shading Factor applies to all tables in a row, regardless of proximity to the edge of the row, as though each row has an infinite number of tables. Proximity to the front of the array is still included, so the correct row count should be used. In order to consider the shading effects of surrounding DC fields and more accurately model shading losses for a DC Field as part of a larger array, this can be applied by activating Array-Based Shading in PlantPredict.

11. Reject all shading factors where the sun is found to be behind the plane of the array, i.e. the back of the module may be lit if the solar incidence angle is less than zero:

The shading factor here is assigned to a value of -1 in order to distinguish it from 0 (fully shaded) and 1 (unshaded). A negative shading factor does not have any effect on the irradiance and will not contribute to the “shading loss” calculation.

The shadows cast by surrounding tables are each determined independently based on their coordinate locations relative to the shaded table. Three cases exist: Transverse shading, which is caused by shading between tables across rows; Adjacent shading, which is caused by shading between tables within rows; and Diagonal shading, which is a combination of transverse and adjacent shading and is caused by tables located diagonally- neither within a row nor directly across rows. The formulas in the algorithm specific to each case are shown in the corresponding sections below identified by the relevant algorithm step.

The net shading fraction on the table is determined by evaluating the total area of the table covered by any of the shadows. This is calculated using an algorithm that determines the area affected by the shadows, accounting for the potential overlapping area of multiple shadows on the shared table.

The shadows considered in the current PlantPredict implementation are due to the tables shown in Figure 22 below, where each of the diagonal tables are in a dark or light red border.

### Figure 22. Tables from which shadows onto the shaded table are considered

Light green and dark green tables are considered depending on the tilt of the tables such that only the tables faced by the shaded table are considered. Light and dark blue tables are considered depending on the slope direction such that only the higher table is considered, since the lower table would never cast a shadow on front of the shaded table. Light and dark red-bordered tables are considered depending on the direction of the sun. These sets of tables are considered based on the sun azimuth angle relative to the azimuth angle faced by the tables.

### Transverse (Row-to-Row) Shading

To calculate the shadow, we need to know the vertical distance between the first row’s top edge and the second row’s top edge. This spacing is defined as the glass-to-glass height difference Δggrr and is a function of the post-to-post spacing Δpp and the slope of the ground in the row-to-row direction 𝜓.

The tilt angle will always be positive for a fixed tilt system due to the user input restrictions. The primary row is always the one faced by the shaded row, which does not change for a fixed tilt system, but flips between rows (also known as columns) on either side of the shaded row, according to the tilt, for a tracker system. Therefore, for a tracker, the shadow of the table in the sun direction is considered based on the tilt by using the negative Δggrr for a negative tilt.

To calculate the shadow, we need to know the vertical distance between the first table’s top edge and the second table’s top edge, which can be calculated between any points along them since they are horizontal. This spacing is defined as the glass-to-glass height difference ΔggTT, and is a function of the table length LT, table-to-table spacing Lg, and the slope of the ground in the table-to-table direction Φ.

Two adjustment factors, AdjFactor and bb1, are used to determine whether the primary table is in the positive or negative y-direction from the shaded table. This logic allows the algorithm to consider only the higher table along the row, which is the only one that can cast a shadow on the shaded table regardless of the sun direction.

For a tracker system,

AdjFactor = -1 and bb1 = 0, for slope Φ ≥ 0

For a fixed tilt system,

AdjFactor = -1 and bb1 = 0, for slope Φ < 0

Otherwise,

AdjFactor = 1 and bb1 = -1

### Diagonal Table Shading

The additional tables considered as casting shadows are located diagonally in the array from the shaded table (neither within a row nor within the column of table). The primary table coordinates Pp are found by a combination of the two coordinates used in Adjacent and Transverse shading. An example diagonal table Pp, in the first row and first column away from the shaded table, is shown below.

It is only necessary to consider diagonal tables in the direction of the sun relative to the shaded table, since only those will cast their shadow onto that table. The algorithm determines whether the primary table is in the positive or negative y-direction from the shaded table, and the corresponding z-direction, by evaluating whether the azimuth of the sun, γs, is greater than the azimuth of the array, γ.

As shown in Step 8 in the general shading algorithm for finding the location of a shadow within the bounds of the shaded table, the location of overlap of two rectangles can be found. This methodology can be used to determine the area over which two shadows overlap each other based on the coordinates of the vertices of each shadow. Since only the total shaded area is important, the overlapping area must not be counted twice as shaded area.

To determine the total area shaded by an arbitrary number of shadows, an algorithm is devised to determine the overlap due to each permutation of the shadows. The total area is calculated as the sum of the overlap area calculated for each permutation made up of an odd number of shadows minus the sum of the overlap area calculated for each permutation made up of an even number of shadows.

# Near Diffuse Shading Factor

Loss of irradiance due to shading is an important calculation for fixed tilt and tracker PV power plants. In addition to modeling direct beam shading, PlantPredict models diffuse shading. Diffuse shading is modeled using a diffuse sky-shading factor, which represents the effective fraction of sky dome lost due to near objects. These calculations are currently restricted to row-to-row shading.

The diffuse shading algorithm used by PlantPredict assumes a clear sky isotropic sky dome. It also assumes that rows of modules have infinite length. Diffuse shading is calculated for the midpoint of the table, represented by Point A in Figure 23, which is deemed representative for the entire table.

## Inputs

*Found using equations from Near Direct Beam Shading

## Algorithm

To estimate near diffuse shading factor, the portion of the sky dome visible to Point A is determined with and without the preceding row of modules. As seen in the equation below, the ratio of the former to the latter is the near diffuse shading factor (i.e. the unshaded fraction).

AreaSky-RTR is the fraction of the sky visible to point A after row-to-row shading. It is not the area of the sky lost due to row-to-row shading. AreaSky-noRTR is the fraction of the sky dome that is visible to point A if the preceding row of modules was not present. Please note that, unless the module tilt = 0°, the denominator of the equation above is not equal to the full sky dome. This is because some of the sky dome is lost to the ground.  The geometry used to define both AreaSky-RTR and AreaSky-noRTR is shown in Figure 24.

### Figure 24. Geometry definitions for diffuse shading model.

In Figure 24, β is equal to the module tilt. γ is the angle between the horizontal running through point A and the line connecting point A to point C (the top of the preceding row of modules). α is the angle between the surface normal from point A and the line connecting point A to point C.

Angles are defined as:

The final equation for AreaSky-RTR is shown as:

Where Atop = tan (β + γ).

Because of the periodicity of the tangent function, if α < 0 then AreaSky-RTR actually becomes the shaded area of the sky, as opposed to the visible portion of the sky, and the following adjustment is made to the output of the equation for AreaSky-RTR.

The final equation for AreaSky-noRTR is nearly identical to that used for AreaSky-RTR and is included below.

Where Aβ = tan (β).

Note that the equations above include the ASHRAE incidence angle modifier b0. PlantPredict uses the same equations to compute diffuse shading and diffuse incidence angle losses. When computing diffuse shading losses, b0 is assumed to be zero. When computing the diffuse incidence angle modifier, μD,IAM, the equation below is used.

The effective sky diffuse irradiance on the tilted plane after diffuse IAM and shading, GD,Shd,IAM, is shown as:

In the above equation, GD,POA is tilted diffuse sky irradiance.

### Sloped Ground Diffuse Shading

For a sloped array, the function is changed such that γ is now a function of the ground slope.

The diffuse shading algorithm accounts for ground slope that would place one row of modules higher, or lower, than the proceeding row. However, it does not currently account for ground slope along the length of rows. Examples of ground slope shading modeled by PlantPredict are shown in Figure 25 through Figure 27.

### Figure 27: Upward sloping ground, when shading is caused by the preceding row.

This algorithm calculates the effective irradiance from the ground to the module corrected for incidence angle losses for a uniform array (tracker or fixed tilt).

## Algorithm

### Single Array

The irradiance from the ground onto a surface tilted with angle β, where GS,POA is the irradiance on the surface, is:

When modeling the ground reflected irradiance on a utility scale solar array, the above equation is insufficient because it does not account for the IAM response of the module surface when receiving ground reflected irradiance. Also, it does not account for row-to-row “ground shading”. The subsequent sections adapt the above equation to account for IAM and ground shading effects.

### Figure 28. Geometry for diffuse ground shading calculation

Consider a two table cross section of a large array as shown in Figure 28. The total irradiance onto the segment between two tables is the sum of a diffuse component and a direct component. The diffuse component is found by calculating the sky view factor between the arrays. The direct component is found using the sun and array angles.  Rewriting the equation for the ground-reflected irradiance for a single plate to include both the diffuse and direct terms results in the following:
Where:
Gs,POA-1D  is the irradiance reflected from the ground to the table surface.
GD,1D is the diffuse irradiance on the segment.
GB,1D  is the beam irradiance on the segment.

### GD, 1D, Diffuse

The total diffuse irradiance on a module in a regular array can be approximated by integrating the product of sky-view factor and IAM factor over a reasonable range in front of the module.

We assume that the reasonable range is from a point on the row in front of the module (-Δpp where Δpp is row-to-row spacing) to the bottom field of view of the module (the closest point on the ground visible to the module), denoted by P, where P = -k cot β.

### Figure 30. Additional angles and segments used in the geometric calculations.

In the figures above, point O is the origin. The coordinate system assumes that right and up are positive directions. Thus all values of x, are negative. Although not shown in the figures, WT is the table width (the distance between points E and B).  Below are the calculations of the geometry of the system:

The coordinates are:

And the angles are:

If λ(x) < 0 then:

The sky-view factor is then:

The fractional sky-view factor is then:

For simplicity, we are using the ASHRAE IAM model:

Where the incidence angle,  is:

### GB,1D, Direct

In this approximation, only the direct beam hitting the ground in front of the module is considered. The direct beam is a step function dependent on the sun angle with respect to the array. Consequently, the limits of integration are a function of solar and array position.

### Reported Albedo Shading

After the separate direct and diffuse albedo shadings are calculated, they are combined into a single reported loss fraction, fL, which is defined as fraction of the single array solution.

# Soiling

Module soiling by dust, snow, or other foreign materials will reflect and scatter some incoming irradiance and render it unavailable for photovoltaic conversion. The simplest approach is to model this as monthly irradiance loss coefficients that affect all irradiance components (diffuse & direct) with equal weight. The soiling level is defined by a vector of loss percentages, one for each month, as defined in the prediction’s Environmental Conditions.

### Figure 31. Example Soiling Loss by Month

These loss factors are then translated into a soiling loss factor, or “soiling ratio”, μdirt,t which is a fraction close to unity, and is the complement of the soiling loss percentage.

For example, a loss coefficient of 0.99 corresponds to a 1% loss of irradiance due to soiling on the modules.

## Algorithm

The soiling factor can be derived by three methods. All inputs indexed by i are defined as monthly values in the Environmental Conditions. The following algorithm assumes a loop through all timesteps, indexed by t.

## Switching the Soiling Model Selection

1.) Case Off — Degenerate case: soiling factor is set to unity for all prediction timesteps regardless of the contents of the site (soiling loss percentage, etc.)

2.) Case On — Base case: soiling factor for the current prediction timestep is based on a lookup of the current monthly soiling loss, as defined in the monthly parameters of the Prediction Environmental Conditions.

# Reflection and Refraction

This routine computes losses due to optical effects at the module level, namely reflection and refraction losses on the sunny-side glass superstrate. This is known as the incidence angle modifier (IAM). Three models are supported: Sandia, ASHRAE, and Tabular.

# Sandia IAM

This model computes the incidence angle modifier coefficients based on the Sandia method.

## Algorithm

Only the beam component ΘEff,D is treated here. Furthermore, it was observed that the best fit of the IAM coefficients for both FSLR anti-reflective coating (ARC) and regular non-ARC modules shows an oscillation below incidence angles of 34°, as shown in Figure 32.

### Figure 32. Comparison of Custom IAM Profile with Sandia Polynomial, Showing Oscillation for Incidence Angles Lower than 34°

To avoid this artifact, the polynomial is clipped to force the IAM to unity, as follows:

## Reference

King, D. L., Kratochvil, J.A., Boyson, W.L., Measuring Solar Spectral and Angle-of-Incidence Effects on Photovoltaic Modules and Solar Irradiance Sensors. Sandia National Laboratories, Albuquerque, NM, September 1997.

# ASHRAE Direct Beam IAM

This model uses the incidence angle modifier coefficients based on the ASHRAE parameterization.

## Algorithm (Direct Beam Component)

1.) Compute the direct beam incidence angle attenuation coefficient as a function of the incidence angle.

2.) The incidence angle function approaches a singularity and is only valid at incidence angles less than 87°.

3.) Reject all shading factors where the sun is found to be behind the plane of the array, i.e. the back of the module may be lit if the solar incidence angle is less than zero:

The shading factor here is assigned to a value of -1 (in order to distinguish it from 0 (fully shaded) and 1 (unshaded). A negative shading factor does not have any effect on the irradiance and will not contribute to the “shading loss” calculation.

# Tabular IAM

This model uses incidence angle modifier coefficients defined in a table in the Module file.

## Algorithm (Direct Beam Component)

1.) Calculate the direct beam incidence angle modifier by doing a cubic spline interpolation on IAM tabular values.

2.) If μIAM,B is greater than 1, set it to 1.

# Diffuse IAM losses

Based on the work of Brandemuehl and Beckman (1980), the following equation integrates the effect of the Angle of Incidence (AOI) for all angles within the field of view to provide an overall AOI correction factor for the diffuse radiation, Fd:

where ω is the solid angle of the incident diffuse irradiance and A is the range of ω. The radiation’s F(AOI) is weighted by its contribution to the in-plane irradiance.

F(AOI) is in this case the tabular IAM values that are part of the module file. The equivalent numerical equation that is being used becomes:

Please see the below drawing for a more detailed explanation:

Reference:

B. Marion “Numerical method for angle-of-incidence correction factors for diffuse radiation incident photovoltaic modules”, Solar Energy, Volume 147, Pages 344-348. 2017. DOI: 10.1016/j.solener.2017.03.027

# Spectral Behavior

The spectral shift factor is representative of how much the performance of a PV system will vary from the nameplate due to differences in outdoor spectrum from flash spectrum (AM 1.5/ASTM G173). In the context of a PV modeling tool, the spectral shift factor can be thought of as a loss factor similar to DC health or soiling. It varies throughout the course of the year and can be positive (a gain instead of a loss).

The spectral shift factor can be calculated in several ways, depending on the availability of weather data and PV module technology. Historically, the air mass portion of Sandia National Lab’s Photovoltaic Array Performance Model has been used to predict the spectral shift for c-Si modules. More recently, several studies have shown that for CdTe modules, this module is highly inaccurate and a model based on precipitable water is instead more appropriate. PlantPredict defaults to the Sandia air mass model for c-Si modules and the precipitable water model for CdTe modules.

The user can override these defaults by using the Spectral Override function in the Environmental Conditions. In the case that this option is selected, PlantPredict will not use the functions below and will use the static monthly values filled in by the user. The user can also turn off spectral calculations in the Advanced Simulation Settings section of the Prediction.

PlantPredict offers a choice of two spectral models in addition to the option to use monthly override or none:

1.) A 1-variable precipitable-water or Sandia model; (PWat)

2.) A 2-variable precipitable-water & air-mass model. (Am_PWat)

The 1-variable spectral response model can be selected at the module-level by a 4-way module-level switch {Off, Sandia, FS Series 4 & Earlier, and FS Series 4-2 & Later} as shown in Figure 33. The Sandia model requires the a0-a4 response coefficients defined under the “advanced parameters” of the module definition. A second set of response coefficients for the Am-PWat model, b0-b5, defines the behavior specifically for the Am-PWat model.

# Sandia Method

This routine computes losses due to spectral response at the module level. This model uses the air mass coefficients computed by Sandia, and is referred to as the f1 (AMa) function in the Array Performance Model. The model is based on an empirically determined polynomial relating the solar spectral influence on the short circuit current to air mass variation over the day.

## Reference

King, D. L., Kratochvil, J.A., Boyson, W.E., Measuring Solar Spectral and Angle-of-Incidence Effects on Photovoltaic Modules and Solar Irradiance Sensors. Sandia National Laboratories, Albuquerque, NM, September 1997.

King, D. L., Boyson, W.E., Kratochvil, J.A., Photovoltaic Array Performance Model. Sandia National Laboratories, Albuquerque, NM, December 2004.

# Precipitable Water Method

The spectral shift factor for CdTe is primarily driven by the amount of precipitable water in the atmosphere, and consequently, the primary input to the calculation is precipitable water. There are several viable sources for precipitable water data: TMY3 data, free online websites, and approximations derived from either relative humidity and ambient temperature data or dewpoint temperature and ambient temperature data.

PlantPredict automatically chooses the best input data source for you, depending on what fields are available in your weather file. If precipitable water data is directly available in your weather file, PlantPredict will use that as input to the calculation. If precipitable water is not available, PlantPredict will estimate precipitable water from relative humidity and ambient temperature. Finally, if no relative humidity data is available, PlantPredict will estimate precipitable water from dewpoint temperature and ambient temperature.

## Algorithm

1.) If neither relative humidity nor precipitable water data is available in the weather file, first estimate the relative humidity from dewpoint and ambient temperatures using the following approximation:

2.) If relative humidity is available in the weather file or if dewpoint was used to estimate relative humidity, estimate the precipitable water using the following approximation, where the air temperature is given in Kelvin, and the relative humidity is given as a percentage:

3.) Find the spectral shift factor:

a.) For First Solar Series 2, Series 3, and Series 4 modules:

b.) For First Solar Series 4-2 and newer modules:

# Precipitable Water & Air Mass Spectral Shift Method

Lee et al, showed that the spectral shift of CdTe and c-Si are both dependent on precipitable water and absolute air mass. In CdTe PV modules, precipitable water is the primary driver of spectral shift, with air mass being of secondary importance. For c-Si PV modules, the opposite is true, and air mass is the primary variable effecting spectral shift.

## Algorithm

Compute the spectral shift factor.

Note that this equation becomes unstable as UW approaches 0 and the air mass gets large near sunrise & sunset, so the author suggests clamping these values as follows:

If monthly spectral shift factors are desired (e.g. to automatically populate monthly factors for a site), compute the irradiance-weighted monthly values from hourly or sub-hourly values. The subscript n denotes the month of the year; subscript i denotes the hourly (or sub-hourly) value within each month.

For example, if the January spectral shift factor is desired, sum the product of all hourly spectral shift factors and global horizontal irradiance for that month, and divide by the sum of the global horizontal irradiance for that month.

## References

Accurate Measurement, Using Natural Sunlight, of Silicon Solar Cells by William M. Keogh and Andrew W. Blakers, http://cses.cecs.anu.edu.au/files/Natural_sunlight_PiPV.pdf

L. Nelson, M. Frichtl, A. Panchula, Changes in Cadmium Telluride Photovoltaic System Performance due to Spectrum by, IEEE Journal of Photovoltaics, Vol. 3, No. 1, January 2013.

M. Lee, L. Ngan, J. Sorenson and A. Panchula, Understanding Next Generation Cadmium Telluride Photovoltaic Performance due to Spectrum, 42nd IEEE Photovoltaic Specialists Conference, June 2015.

Lee, M, Panchula, A., Spectral Correction for Photovoltaic Module Performance Based on Air Mass and Precipitable Water, 43rd IEEE Photovoltaic Specialists Conference, June 2016