Solar irradiation from the energy production of residential PV systems

Accepted Manuscript Solar irradiation from the energy production of residential PV systems Cédric Bertrand, Caroline Housmans, Jonathan Leloux, Michel Journée PII:

S0960-1481(18)30181-2

DOI:

10.1016/j.renene.2018.02.036

Reference:

RENE 9767

To appear in:

Renewable Energy

Received Date: 24 October 2017 Revised Date:

11 January 2018

Accepted Date: 6 February 2018

Please cite this article as: Bertrand Cé, Housmans C, Leloux J, Journée M, Solar irradiation from the energy production of residential PV systems, Renewable Energy (2018), doi: 10.1016/ j.renene.2018.02.036. This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT

RI PT

Solar irradiation from the energy production of residential PV systems C´edric Bertrand †∗, Caroline Housmans †, Jonathan Leloux ‡, and Michel Journ´ee † † Royal Meteorological Institute of Belgium, Brussels, Belgium.

SC

‡ Universidad Polit´ ecnica de Madrid, Madrid, Spain

M AN U

Abstract

Considering the dense network of residential photovoltaic (PV) systems implemented in Belgium, the paper evaluates the opportunity of deriving global horizontal solar irradiation data from the electrical energy production registered at PV systems. The study is based on one year (i.e. 2014) of hourly PV power output collected at a representative sample of roughly 1500 residential PV installations. Validation is based on ground-based measurements of solar radiation

TE D

performed within the network of radiometric stations operated by the Royal Meteorological Institute of Belgium and the method’s performance is compared to the satellite-based retrieval approach.

EP

Our results indicate that the accuracy of the derived solar irradiation data depends on a number of factors including the efficiency of the PV system, the weather conditions, the density of PV systems that can be used for the tilt to

AC C

horizontal conversion, other data sources that can be accessed to complement the PV data. In particular, the computed solar irradiation data degrade as the information about the orientation and tilt angles of the PV generator becomes more inaccurate. It is also found that there are certain sun positions (i.e. low solar elevations) for which the method fails to produce a valid estimation. ∗ Corresponding

author. Tel.:+32 2 373 05 70 Email address: [email protected] (C´ edric Bertrand † )

Preprint submitted to Renewable Energy

January 8, 2018

ACCEPTED MANUSCRIPT

Keywords: Photovoltaic system, solar radiation, decomposition model,

1

RI PT

transposition model, remote sensing, ground measurements, crowdsourcing

1. Introduction

Appropriate information on solar resources is very important for a variety

3

of technological areas, such as: agriculture, meteorology, forestry engineering,

4

water resources and in particular in the designing and sizing of solar energy

5

systems. Traditionally, solar radiation is observed by means of networks of me-

6

teorological stations. However, costs for installation and maintenance of such

7

networks are very high and national networks comprise only a few stations.

8

Consequently the availability of solar radiation measurements has proven to be

9

spatially and temporally inadequate for many applications.

M AN U

SC

2

10

Over the last decades, satellite-based retrieval of solar radiation at ground

12

level has proven to be valuable for delivering a global coverage of the global

13

solar irradiance distribution at the Earth’s surface (e.g. [1], [2] [3], [4]). The

14

recent deployment of solar photovoltaic (PV) systems offers an opportunity to

15

extend the in-situ characterization of the incoming solar radiation at the Earth’s

16

surface through the reverse engineering conversion of the energy production of

17

the PV systems to global solar irradiation data.

EP

18

TE D

11

During the last years, thanks to the newest developments in telecommunica-

20

tion technologies, an increasing amount of electrical (power, frequency, current,

AC C

19

21

voltage etc ...) data from the operation of PV systems has has become avail-

22

able. These operational data are now typically available at time intervals of 10

23

minutes, opening the doors to a totally new set of data analysis procedures that

24

were difficult to imagine just a few years ago. Smart energy meters are pro-

25

jected to become widespread in a couple of years, which will enable to monitor

26

the energy production and other electrical outputs of millions of PV systems

27

in the world. Therefore, in the near future, it will be possible to know the

2

ACCEPTED MANUSCRIPT

energy production of PV systems located at millions of points spread over the

29

whole globe. This conveys the potential that solar irradiation data can also be

30

retrieved at all these points.

RI PT

28

31

However, using the energy production registered at PV systems as a solar

33

irradiation sensor is not straightforward (e.g. [5], [6]). It requires first to derive

34

the solar irradiation from the energy production of the PV system (knowing

35

that the power output of a PV system is not directly proportional to the solar

36

irradiance that it receives). Second because modules are installed at a tilt angle

37

close to local latitude to maximize array output (or at some minimum tilt to

38

ensure self-cleaning by rain) this requires to convert the retrieved tilted global

39

solar irradiance to horizontal. Towards this objective, operational data from a

40

representative sample of Belgian residential PV installations (note that tilt and

41

azimuth of residential PV systems often depend on the particular rooftop on

42

which they are installed, rather than being designed for optimal performance )

43

have been considered to assess the performance of such an approach.

M AN U

TE D

44

SC

32

In Belgium, the total installed PV capacity has increased dramatically in

46

recent years from 102.6 MW in 2008 (26.55 MW in 2007) to 3423 MW at the

47

end of 2016, according to data collected by local renewable energy association

48

APERe (http://www.apere.org/), which has combined the figures released by

49

the country’s three energy regulators Brugel (https://www.brugel.brussels/),

50

VREG (http://www.vreg.be/)and CWaPE (http://www.cwape.be/). By con-

51

trast to the multi-megawatt PV plants installed in the South of Europe, res-

AC C

EP

45

52

idential and commercial rooftop PV systems with less than 10 kW capacity

53

represent the majority of the installed capacity in Belgium. As an illustration,

54

in 2016, the country installed about 170 MW across 25.000 PV systems.

3

55

2. Data

56

2.1. Residential PV systems data

RI PT

ACCEPTED MANUSCRIPT

This work is based on one year (i.e. 2014) of hourly PV power output

58

collected at more than 2893 PV systems in Belgium installed from 2008 to

59

2013. PV generation data was collected via the Rtone company website (Rtone,

60

http://www.rtone.fr). The PV energy production data provided by Rtone was

61

monitored using the commercial Rbee Solar monitoring product, which measures

62

the energy production with a smart energy meter at a 10-min time interval. The

63

information concerning the PV systems (i.e. metadata) were supplied by their

64

owners. Each PV system is localized by its latitude and longitude, completed

65

with the corresponding altitude. The PV generator is characterized by its ori-

66

entation and tilt angles, its total surface, and its total peak power. Additional

67

information about the PV module manufacturer and model, inverter manufac-

68

turer and model, installer, year of installation, PV cell/module technology can

69

also be provided.

TE D

70

M AN U

SC

57

71

Metadata has been subjected to several checks in order to isolate and remove

72

as much erroneous information as possible. The standard set of filters employed

73

before analyses are the following:

75

1. selection of single array systems since generation data cannot be decom-

EP

74

posed into constituent arrays, 2. selection of systems located in Belgium,

77

3. selection of systems with and orientation between -90o and +90o from

AC C

76

78

south and a tilt from horizontal smaller than 60o .

79

In some cases, system details were investigated manually to verify to a good

80

degree of confidence whether the systems orientation, tilt or peak power are

81

correct. Indeed, metadata are prone to errors on part of the owner, by entering

82

incorrect system parameters but also in some cases due to installer’s conventions.

83

As an example, tilt angle can be erroneously reported as 0o simply because 0 is

4

ACCEPTED MANUSCRIPT

used as default value in case of missing information. Another identified limita-

85

tion is that many angles values are rounded as multiples of 5o .

RI PT

84

86

1470 systems have been selected and considered in this study after these

88

preliminary checks (see Figure 1 for the location of the selected PV systems).

89

It is worth pointing out that due to availability reasons, most of the data comes

90

from Wallonia and Brussels. As indicated in Table 1, the vast majority of the

91

selected PV generators have a tilt angle between 20o and 50o , which generally

92

corresponds to the slope of the roofs on which they are installed ([7]).

SC

87

94

M AN U

93

2.2. Ground stations measurements

The Royal Meteorological Institute of Belgium (RMI) has a long term ex-

96

perience with ground-based measurements of solar radiation in Belgium (unin-

97

terrupted 30 min average measurements in Uccle since 1951, in Oostende since

98

1958, and in Saint-Hubert since 1959). Uccle is one of the 22 Regional Radiation

99

Centres established within the WMO Regions. The incoming global horizontal

100

solar irradiance is currently measured in 14 Automatic Weather Stations (AWS)

101

in addition to the measurements performed in our main/reference station in Uc-

102

cle (see Table 2).

103

TE D

95

At the reference station in Uccle, measurements of the global horizontal so-

105

lar irradiance (GHI) are performed by Kipp & Zonen CM11-secondary standard

106

pyranometers. The CM11 pyranometer has a thermopile detector and presents

AC C

EP

104

107

a reduced thermal offset (i.e. about 2 Wm−2 at 5 K/h temperature change).

108

The pyranometer directional error (up to 80o with 1000 Wm−2 beam) is less

109

than 10 Wm−2 and its spectral selectivity (300-1500 nm) is smaller than 2%.

110

At the 14 RMI’s AWS, GHI observations for the year 2014, were made with a

111

Kipp & Zonen CNR1 net radiometer. It consists of a pyranometer (model type

112

CM3 complying with the ISO Second Class Specification) and a pyrgeometer

113

(model type CG3) pair that faces upward and a complementary pair that faces 5

ACCEPTED MANUSCRIPT

downward. The pyranometers and pyrgeometers measure shortwave and far in-

115

frared radiation, respectively. The CM3 pyranometer has a thermopile sensor

116

and presents a limited thermal offset (i.e. ± 4 Wm−2 at 5 K/h temperature

117

change). Its directional error (with 1000 Wm−2 beam) is ± 25 Wm−2 and its

118

spectral selectivity (350-1500 nm) is ± 5%.

119

RI PT

114

All ground measurements are made with a 5-s time step and then integrated

121

to bring them to a 10-min time step. The 10-min data have undergone a series

122

of automated quality control procedures ([8]) prior to be visually inspected and

123

scrutinized in depth by a human operator for more subtle errors. Because our

124

radiometric station in Saint-Hubert (AWS 6476) was known to operate defi-

125

ciently in 2014, all GHI measurements from this station were discarded. The

126

geographical location of the remaining 13 ground measurement sites is provided

127

in Figure 1 together with the selected residential PV systems. Note that be-

128

cause the data quality control revealed that the CNR1 net radiometer installed

129

in the Buzenol station (i.e. AWS 6484) has only performed well intermittently

130

during the year 2014, GHI measurements from this station were not used for

131

validation purpose.

132

TE D

M AN U

SC

120

3. Methods

134

3.1. Conversion of PV system energy production to tilted global solar irradiation

135

The initial step of the approach consists in the derivation of global irradi-

AC C

EP

133

136

ance in plane of array, Gt , from the specific power output, P , of a PV system.

137

Numerous models to calculate P from Gt exist in the literature (e.g. [9], [10]).

138

However, it is well known that the energy conversion efficiency of PV mod-

139

ules depends on a number of different influences. Losses in PV systems can be

140

separated in capture losses and system losses (e.g. [11], [12]). Capture losses

141

are caused, e.g., by attenuation of the incoming light, temperature dependence,

142

electrical mismatching, parasitic resistances in PV modules and imperfect max-

6

ACCEPTED MANUSCRIPT

imum power point tracking. System losses are caused, e.g., by wiring, inverter,

144

and transformer conversion losses. All these effects cause the module efficiency

145

to deviate from the efficiency measured under Standard Test Condition (STC),

146

which defines the rated or nominal power of a given module.

147

149

According to [13], the direct current (DC) power output of a PV generator can be properly described by:

SC

148

RI PT

143

150

PDC = P ?

    Gef f Gef f Gef f ? . 1 + κ (T − T ) . a + b + c ln . fDC (1) c c G? G? G? ?

refers to STC [Irradiance: 1000 W m−2 ; Spectrum: AM

where the symbol

152

1.5; and cell temperature: 25o C], PDC is the DC power output of the PV gen-

153

erator (W), P ? is the nameplate DC power of the PV generator (i.e. power at

154

maximum-power point, in W), Gef f is the effective global solar irradiance (W

155

m−2 ) received by the PV generator (it takes into consideration the optical ef-

156

fects related to the solar incidence angle), G? is the global solar irradiance under

157

STC (W m−2 ), κ is the coefficient of power variation due to cell temperature

158

(%/o C) , Tc and Tc? are respectively the cell temperatures under operating and

159

STC conditions (o C), the three parameters a, b and c describe the efficiency

160

dependence on irradiance, and fDC is a coefficient that lumps together all the

161

additional system losses in DC (e.g. technology-related issues, soiling and shad-

162

ing).

TE D

EP

163

M AN U

151

The first round brackets on the right-hand side of Eq. 1 takes into account

165

the dependence of the efficiency with the temperature. The formulation goes a

166

long way back ([14], [15]) and considers that the PV module efficiency is affected

167

by temperature, decreasing at a constant rate. Handling with this model for-

168

mulation just requires standard information: P ? is the PV array rated power,

169

which can be estimated as the product of the number of PV modules constitut-

170

ing the PV array multiplied by their nameplate STC power, and κ is routinely

171

measured in the context of worldwide extended accreditation procedures: IEC

AC C

164

7

ACCEPTED MANUSCRIPT

Standard-61215 (2005) and IEC Standard-61646 (2008) for crystalline silicon

173

and thin film devices, respectively. P ? and κ values are always included in PV

174

manufacturer’s data sheets or in more specific information as flash-reports.

RI PT

172

175

The second round bracket on the right-hand side of the Eq. 1 considers the

177

efficiency dependence on irradiance. That was initially attempted by adding a

178

base 10 logarithm ([14]) but it is better implemented by this empirical model

179

proposed by [16] and [17] where a, b and c are empirical parameters. The effi-

180

ciency increases with decreasing irradiance, due to series resistance effects, are

181

represented by the term (b . Gef f /G? ), providing b ≤ 0, while the efficiency

182

decreases with decreasing irradiance, due to parallel resistance effects, are rep-

183

resented by the term (c ln Gef f /G? ), providing c ≥ 0.

M AN U

SC

176

184

185

186

The corresponding alternating current (AC) power output of the PV system from this DC power at the inverter entry is given by:

187

TE D

PAC = PDC ηIN V fAC

(2)

where PAC is the AC power output of the PV generator, ηIN V is the yield of the

189

inverter, and fAC is a coefficient that lumps together all the technology-related

190

additional AC system losses.

191

192

The energy produced during a period of time T is finally given by:

AC C

193

EP

188

Z

t=T

EAC =

PAC dt

(3)

t=0

194

To assess the technical quality of a particular PV system, energy performance

195

indicators are obtained by comparing its actual production along a certain pe-

196

riod of time with the production of a hypothetical reference system (of the same

197

nominal power, installed at the same location, and oriented the same way). The

198

Performance Ration (PR) which is the quotient of alternating current yield

199

and the nominal yield of the generators direct current, is by far the most widely

8

ACCEPTED MANUSCRIPT

200

used performance indicator. It is defined mathematically as:

PR =

ηachieved EAC /Gt = ? ? ηspec PN /G

RI PT

201

(4)

where PN? is the nominal (or peak) DC power of the PV generator, understood

203

as the product of the number of PV modules multiplied by the corresponding

204

in-plane STC power. Because EAC , PN? and Gt are given by the billing energy

205

meter of the PV installation, the PV manufacturer and the integration of a

206

solar irradiance signal, the PR value can be directly calculated. The difference

207

between 1 and PR lumps together all imaginable energy losses (i.e. capture

208

losses and system losses).

M AN U

209

SC

202

For a given PV system and site, the P R value tends to be constant along

211

the years, as much as the climatic conditions tend to repeat. When sub-year

212

periods are considered, the P R dependence on unavoidable and time-dependent

213

losses requires corresponding correction in order to properly qualify the tech-

214

nical quality of a PV system. Based on Eqs. 1 to 3, we can reformulate Eq. 4 as:

215

TE D

210

1 .PR = 1 fG . fT . fAC . fP DC . fBOS

216

(5)

in which the losses have been lumped into five main categories: 1. fG : PV module’s yield in function of incident irradiance level,

218

2. fT : PV module’s yield in function of cell’s temperature,

219

3. fAC : yield of the conversion from DC to AC.

220

4. fP DC : yield that represents the ratio of the real DC power and the rated

AC C

EP

217

221

222

223

224

DC power,

5. fBOS : yield of the balance of system. Three of these five losses parameters can be expressed analytically. Based

on Eq. 1, the efficiency dependence on irradiance is:

225

fG = a + b

Gef f Gef f + c ln ? G G? 9

(6)

ACCEPTED MANUSCRIPT

However, such a formulation of the fG parameter is useless here since the

227

effective irradiance, Gef f , is by definition unknown in our case. To overcome

228

such a limitation, fG is split into its two main contributing factors:

RI PT

226

229

fG = firr . finc

(7)

where firr represents the variation in the PV module efficiency with the level

231

of the solar irradiance and finc the variation in the PV module efficiency as a

232

function of the incidence angle of the solar irradiance, respectively. Then, ap-

233

proximating the ratio Gt /G? by the Capacity Utilization Factor (CUF) defined

234

as:

M AN U

SC

230

235

CU F =

EAC T . PN?

(8)

236

with EAC the energy produced during a period of time T (see Eq. 3) and PN?

237

the nominal (or peak) DC power of the PV system, firr can be estimated by:

238

(9)

TE D

firr = a + b . CU F + c ln(CU F )

In this equation, the three parameters a, b and c vary according to the con-

240

sidered PV module technology. Values representative of crystalline silicon cells

241

technology (i.e., a=1, b=-0.01 and c=0.025) have been assumed for all PV mod-

242

ules here. Finally, based on [18], finc can be expressed as follows:

243

EP

239

finc = 1 −

1 − exp(−(cos θi )/αr ) 1 − exp(−1/αr )

(10)

where θi is the irradiance angle of incidence and αr the angular loss coefficient,

245

an empirical dimensionless parameter dependent on the PV module technology

246

and the dirtiness level of the PV module. Typical αr values range from 0.16 to

247

0.17 for commercial clean crystalline and amorphous silicon modules. In this

248

work a value of 0.20 has been assumed for αr which is a typical value for crys-

249

talline silicon PV modules presenting a moderate level of dirtiness.

AC C

244

250

10

ACCEPTED MANUSCRIPT

251

The second factor, fT , is defined as:

RI PT

252

fT = 1 + κ (Tc − Tc? )

(11)

where the operating temperature of the solar cell, Tc , is calculated from the

254

ambient temperature, Ta , using the following equation based on the Nominal

255

Operation Cell Temperature (NOCT) defined as the temperature reached by

256

the cells when the PV module is exposed to a solar irradiance of 800 W.m−2 ,

257

an ambient temperature of 20o C, and a wind speed of 1 m.s−1 (it is obtained

258

from the manufacturer datasheets):

M AN U

259

Tc = Ta + 260

261

N OCT − 20 . Gt 800

(12)

Similarly to Eq. 6 the CUF approximation is used to estimate Tc reformulating Eq. 12 as follows:

262

 Tc = Ta +



(N OCT − 20)/800 . 1000 . CU F

(13)

The third factor, fAC , is computed from:

264

TE D

263

SC

253

fAC =

PAC . ηEU R

PN?

(14)

where, the so-called ”European efficiency” of the inverter, ηEU R , is given by the

266

formula:

267

EP

265

(15)

AC C

ηEU R = 0.03 η5 + 0.06 η10 + 0.13 η20 + 0.1 η30 + 0.48 η50 + 0.2 η100 268

with η5 , η10 , η20 , η30 , η50 , η100 the instantaneous power efficiency values at 5%,

269

10%, 20%, 30%, 50% and 100% load.

270

271

The fourth factor, fP DC , as well as the fBOS factor cannot not be directly

272

estimated because the real energetic behavior of each PV system is a priori

273

unknown. Lumping both factors together into a new losses factor, fP ERF , it

11

ACCEPTED MANUSCRIPT

274

follows from Eqs. 4 and 5 that

fP ERF =

1 fG . fT . fAC

.

EAC /Gt PN? /G?

RI PT

275

(16)

where the fP ERF factor sums up all the performance losses that, on the first

277

hand could be avoided and, on the other hand that cannot be modeled through

278

a simple and general analytical expression. This factor can be estimated for

279

each PV system from historical data of EAC and Gt using the EAC data di-

280

rectly from the energy meters and Gt data obtained from the combination of

281

clear-sky radiative model simulations and cloud cover information. It is worth

282

pointing out that to reduce the uncertainties in its estimation, the fP ERF factor

283

was determined on a monthly basis from clear sky situations.

M AN U

SC

276

284

Similarly to [19], clear-sky situations were determined from the PV systems

286

energy production time series using a modified version of the algorithm devel-

287

oped by [20]. For each PV system, fP ERF , was calculated as being the ratio

288

between the electrical energy produced by the PV system corrected by the three

289

other losses factors (i.e. fG , fT , and fAC ) together with the quotient PN? /G? and

290

the calculated in-plane clear-sky irradiation received by the PV system during

291

the considered month. Clear sky simulations were carried out by running the In-

292

eichen and Perez clear-sky model [21] using monthly mean climatological Linke

293

turbidity values from PVGIS/CMSAF (http://re.jrc.ec.europa.eu/pvgis/). Sim-

294

ulated clear sky global horizontal irradiances were then transposed to tilted

295

clear sky global irradiance using the ERB decomposition model ([22]; see Ap-

296

pendices A.1) and the HAY transposition model ([23]; see Appendices B.1).

AC C

EP

TE D

285

297

298

Finally, once all losses factors are estimated, the derivation of the in-plane

299

hourly global solar irradiance from the hourly PV system energy production is

300

given by:

301

Gt =

1 EAC . fG .fT .fAC .fP ERF PN? /G? 12

(17)

ACCEPTED MANUSCRIPT

where all losses factors except fP ERF are evaluated on a hourly basis using air

303

temperature measurements performed within the RMI’s AWS spatially interpo-

304

lated at each of the PV system locations in the computation of the operating

305

solar cells temperature, Tc (see Eq. 13). fP ERF are determined monthly from

306

the previous month EAC data.

307

3.2. Tilt to horizontal global solar irradiance transposition

RI PT

302

The next step consists in the conversion of the retrieved in-plane global solar

309

irradiance values from the PV systems energy outputs to global horizontal solar

310

irradiance at each of the PV systems location. If many transposition models

311

have been proposed in the literature (see [24] for a review) to convert solar ir-

312

radiance on the horizontal plane, Gh , to that on a tilted plane, Gt , the inverse

313

process (i.e. converting from tilted to horizontal) is only poorly discussed in

314

literature (e.g. [25], [26], [27], [28], [5], [6], [29]). The difficulty relies on the fact

315

that the procedure is analytically not invertible.

M AN U

SC

308

316

318

Transposition models have the general form:

TE D

317

Gt = Bt + Dt + Dg

(18)

where the tilted global solar irradiance, Gt , is expressed as the sum of the in-

320

plane direct irradiance, Bt , in-plane diffuse irradiance, Dt , and the irradiance

321

due to the ground reflection, Dg . The direct component, Bt , is obtained from:

322

EP

319

323

cos θi = Bh rb (19) cos θz with Bn the direct normal irradiance and Bh the direct irradiance on a horizon-

324

tal surface, respectively. θi is the incidence angle and θz the solar zenith angle,

325

respectively. Parameter rb = cos θi / cos θz is a factor that accounts for the di-

326

rection of the beam radiation. The diffuse component, Dt , and the irradiance

327

due to the ground reflection, Dg , can be modeled as follows:

AC C

Bt = Bn cos θi = Bh

328

Dt = Dh Rd 13

(20)

ACCEPTED MANUSCRIPT

Dg = ρGh Rr

(21)

where Dh is the diffuse horizontal irradiance, Gh the global horizontal irradi-

330

ance (i.e., Gh = Dh + Bh ), Rd the diffuse transposition factor and ρ the ground

331

albedo. The transposition factor for ground reflection, Rr , can be modeled un-

332

der the isotropic assumption (e.g. [30]) as follows:

Rriso =

1 − cos β 2

(22)

SC

333

RI PT

329

where, β, is the tilt angle of the inclined surface. See Figure 2 for angles defini-

335

tion.

M AN U

334

336

337

Considering the effective global horizontal transmittance, Kt , the direct

338

normal transmittance, Kn , and the diffuse horizontal transmittance, Kd (i.e.,

339

Kd + Kn = Kt ):

340

Gh = Kt Io cos θz

(23)

TE D

Bn = Kn Io

Dh = Kd Io cos θz

341

where Io is the extraterrestrial normal incident irradiance, Eq. 18 can be rewrit-

342

ten as:

343

EP

     Kd Kd + cos θz Rd + ρRr Gt = Kt Io cos θi 1 − Kt Kt

(24)

It comes from Eq. 24 that when only one tilted global solar irradiance mea-

345

surement is considered, the conversion of Gt to Gh requires the use of a decompo-

AC C

344

346

sition model (i.e. model that separate direct and diffuse solar components from

347

the global one) to estimate Kd from Kt in addition to a transposition model

348

to solve the inverse transposition problem. Eq. 24 is solved by an iteration

349

procedure, varying the target quantity Gh (through Kt ) until the resulting G0t

350

matches the input Gt (e.g. [26], [28], [5], [29]). Note that an alternative method

351

to Eq. 24 based on the Olmo model ([31]) that presents the particularity of be-

352

ing analytically invertible was proposed by [5]. But, if the overall performance 14

ACCEPTED MANUSCRIPT

of the inverted Olmo model was found comparable with the other approach,

354

the results were slightly worse than those obtained by inverting the decompo-

355

sition and transposition models in combination with an iterative solving process.

RI PT

353

356

When two (or more) tilted irradiances values (with different tilt angles

358

and/or orientations) are involved in the inverse transposition process, only a

359

transposition model is required. The idea that simultaneous readings of a multi-

360

pyranometers system can be used to disangle the various components of solar

361

radiation on inclined surfaces was originally proposed by [25] to solve in remote

362

locations the periodic adjustment required by normal incidence and shadow-

363

band pyranometers to ensure that their readings remain accurate when long-

364

term data acquisition is in progress.

M AN U

SC

357

365

Given n tilted pyranometers (with different inclinations and/or orientations),

367

the inverse transposition problem can be represented in the matrix form (e.g.

368

[27]):

369

TE D

366

xT Λ x + B − C = 0

(25)

where Λ = {Ai } is a 2 × n × 2 third-order tensor, B = {Bi } is a n × 2 matrix,

371

C is a column vector with n given entries, and x a column vector with 2 variables:

AC C

EP

370

15

ACCEPTED MANUSCRIPT

372



0

Ai



RI PT

 ∈ <2×n×2 Ai =  Ai 0   C1 B1      C2 B2  n×2  B= ..  ∈ <  ..  . .    Cn Bn   Gt1      Gt2  n  C= .  ∈<  ..   

M AN U

SC

(26)

Gtn



xT = Dh

 2 Bh ∈ <

373

where the coefficients Ai , Bi and Ci depend on the considered transposition

374

model.

375

The least square (hereafter referred to as LS) solution to Eq. 26 is given by:   1 min P (x) = ||xT Λ x + B − C||2 : x ∈ <2 (27) 2

TE D

376

with ||.|| referring to the Euclidean norm. However, the LS is hard to solve and

378

a standard technique to resolve Eq. 27 is to use a Newton type iteration method

379

(e.g. [32]). As an alternative, Eq. 26 can also be solved by minimizing the errors

380

(this approach is hereafter denoted to as EM - Errors Minimization). In this

381

case, the solution is to minimize

AC C

EP

377

382

( min E(x) =

n X

) 2i (x)

2

:x∈<

(28)

i=1

383

where, i (x) = (Ai Dh Bh + Bi Bh + Ci Dh ) − Gti , with i = 1, ..., n denoting the

384

tilted pyranometer.

16

ACCEPTED MANUSCRIPT

385

4. Results The Mean Bias Error (MBE) and the Root Mean Square Error (RMSE)

387

statistical error indexes have been used to evaluate the prediction of the global

388

horizontal solar irradiance from the energy production of residential PV systems. n

1X (ei ) n i=1 v u n u1 X RM SE = t (e2 ) n i=1 i M BE =

390

SC

389

RI PT

386

where ei = (Gi,e − Gi,o ) is the residual value, Gi,e are the estimated values and

392

Gi,o represent the observed measurements. A positive MBE (resp. a negative

393

MBE) means that the model tends to overestimate (resp. underestimate) the

394

observed measurements.

395

M AN U

391

To obtain dimensionless statistical indicators we express MBE and RMSE

397

as fractions of mean solar global irradiance during the respective time interval:

398

TE D

396

M BE ¯ M RM SE RM SE[%] = ¯ M M BE[%] =

399

401

402

¯ = where M

1 n

n P

(Gi,o ) is the measurements mean.

EP

400

i=1

Statistical error indexes were computed between in situ hourly irradiance measurements and the estimations computed from the hourly energy produc-

404

tions of residential PV systems surrounding the measurement stations. An ini-

405

tial radius of 5 km centered on the station was considered to select the residential

406

PV systems for the validation purpose. When less than 4 PV installations were

407

found within the delimited area, the radius was extended to 10 km. Table 3

408

indicates for each of our measurement sites the number of neighboring PV in-

409

stallations used for validation. No PV system was found in the vicinity of the

AC C 403

17

ACCEPTED MANUSCRIPT

Sint-Katelijn-Waver, Retie and Mont-Rigi stations (i.e. AWS 6439, AWS 6464

411

and AWS 6494, respectively) and 3 others stations only have one surrounding

412

residential PV system. At the opposite, the maximum number of installations

413

surrounding a station is 37 for Ernage (i.e. AWS 6459).

414

RI PT

410

Based on our former evaluation of the inverse transposition problem ([29]),

416

two different approaches have been considered to compute the global horizontal

417

solar irradiance from the PV systems energy production. In the first approach,

418

the tilt to horizontal conversion is performed independently at each PV installa-

419

tions surrounding the validation site using Eq. 24 with the OLS decomposition

420

model (see Appendices A.2) and the SKA transposition model (see Appen-

421

dices B.2). The median value of the individual PV system estimates is consid-

422

ered in the validation against AWS observations. This approach is referred to

423

as 1− PV-M hereafter. In the second approach all individual tilted global solar

424

irradiance estimates are used simultaneously and the tilt to horizontal conver-

425

sion is solved by EM (see Eq. 28) using the Powell’s quadratically convergent

426

method ([33]) and the SKA transposition model (see Appendices B.2). Indeed

427

the minimization carried out by using the Powell’s method has been found to

428

systematically outperform the LS solution ([29]). It is a generic minimization

429

method that allows to minimize a quadratic function of several variables with-

430

out calculating derivatives. The key advantage of not requiring explicit solution

431

of derivatives is the very fast execution time of the Powell method. In order to

432

avoid the problem of linear dependence in the Powell’s algorithm, we adopted

433

the modified Powell’s method given in [34] and implemented in [35]. This second

AC C

EP

TE D

M AN U

SC

415

434

approach is hereafter denoted to as X− PV-EM.

435

436

Performance of the two approaches in the GHI estimation from PV systems

437

AC power output was evaluated against in-situ observations and furthermore

438

compared to the performance of GHI estimates retrieved from Meteosat Sec-

439

ond Generation (MSG, [36]) satellite images as implemented on an operational

440

basis at RMI. Description of the RMI’s MAGIC/Heliosat-2 algorithm used to 18

ACCEPTED MANUSCRIPT

retrieve the solar surface irradiance at the SEVIRI imager spatial sampling dis-

442

tance above Belgium (e.g. about 6 km in the north–south direction and 3.3 km

443

in the east–west direction) from MSG images is provided in Appendices C.

RI PT

441

444

While the MSG based retrieval method always provides GHI estimates dur-

446

ing day time, there are certain sun positions for which the PV systems power

447

output method fails to produce a valid estimation. Unsuccessful tilt to horizon-

448

tal conversions are systematically found at solar elevations lower or equal to 5o

449

for both the 1− PV-M and the X− PV-EM approaches irrespective of the number

450

of PV systems involved in the conversion process and their angular configura-

451

tions (i.e. tilt and azimuth angles). Failure rates reported for the 1− PV-M

452

and the X− PV-EM approaches at each validation sites are provided in Table 3

453

together with the total number of available hourly data points at each location

454

for the year 2014. Unsurprisingly the largest failure rates (up to nearly 40% in

455

the case of the Melle station -AWS 6434-) are found at validation sites where

456

only one PV installation is available. With more PV systems, the number of

457

unsuccessful conversions after sunrise and before sunset is decreased. Table 3

458

tends to indicate that 1− PV-M starts to produce valid results at lower solar ele-

459

vation conditions than X− PV-EM (i.e. an overall failure rate of about 12.4% is

460

reported for 1− PV-M and of 19.6% for X− PV-EM, respectively) but it is largely

461

relying on the angular configurations of the PV installations found within the

462

group of PV systems.

M AN U

TE D

EP

AC C

463

SC

445

464

4.1. Hourly validation

465

Table 4 compares hourly GHI estimates derived from PV production data

466

with the 1− PV-M and X− PV-EM approaches as well as from MSG images with

467

the corresponding ground measurements. To ensure that the comparisons are

468

made between comparable data, special attention was given to the coherence of

469

the data, the precision of the time acquisition, and the synchronization of the

470

different data sets with the ground measurements. Because of inaccuracies in the 19

ACCEPTED MANUSCRIPT

orientations and/or inclinations of the PV installations provided by the PV sys-

472

tems installers or owners, GHI computation from the energy production of only

473

one installation can generate RMSE values as large as 189.34 W.m−2 or 57.8%

474

(i.e. at the Middelkerke validation site, AWS 6407). Increasing the number of

475

PV installations involved in the estimation process smoothes the GHI estima-

476

tion to some extent. This is particularly apparent for 1− PV-M which globally

477

presents lower RMSE values than found for X− PV-EM (i.e. an overall RMSE

478

value of 113.5 W.m−2 or 41.4% is reported for 1− PV-M and of 121.9 W.m−2 or

479

44.4% for X− PV-EM, respectively). However, sensitivity experiments in which

480

the number of PV installations involved in the GHI determination was varying

481

revealed a larger variability in the resulting GHI estimations for 1− PV-M than

482

found for X− PV-EM which produces a more stable solution.

M AN U

SC

RI PT

471

483

[5] reported a somewhat similar mean RMSE error of about 40% for GHI es-

485

timates derived from 5 years (from 2010 through 2014) of five-minute resolution

486

records of specific power of 45 PV systems in the region of Freiburg, Germany.

487

In contrast, GHI retrieval from MSG images shows a better performance with

488

an overall RMSE of 55.8 W.m−2 or 20.3%. Moreover, while the satellite re-

489

trieval tends to slightly overestimate the GHI values (i.e. MBE values ranging

490

from 0.16 to 5.87%), there is no clear trend in the GHI computation from mea-

491

sured AC PV output power. Positive and negative biases are reported for both

492

1− PV-M and X− PV-EM approaches. Moreover it can appear that the sign of

493

the bias even differs from one approach to the other (e.g. a negative MBE value

494

of -9.08 W.m−2 or -3.37% is reported at the Humain validation site (AWS 6472)

AC C

EP

TE D

484

495

for 1− PV-M approach while an overestimation of 7.86 W.m−2 or 2.92% is found

496

for X− PV-EM. In general, the magnitude of the bias is lower with X− PV-EM

497

than with 1− PV-M (i.e. MBE values ranging from -5.6% to 5.5% and from

498

-9.9% to 7.8%, respectively).

499

Table 5 compares the performance of the 1− PV-M and X− PV-EM approaches together with the MSG-based retrieval method in the hourly GHI 20

ACCEPTED MANUSCRIPT

computation as a function of the sky conditions. Sky-type classification is made

Kt

Kt0 =

 1.031 exp

500

where Kt =

Gh E

−1.4 0.9+9.4/M

 + 0.1

RI PT

upon the modified clearness index proposed by [37] and defined as

is the usual clearness index (i.e. the horizontal surface so-

503

air mass and γ is the solar elevation angle in degrees (i.e. γ = 90 − θz ). This

504

zenith angle-independent expression of the clearness index, based on Kasten’s

505

([38]) pyrheliometric formula, ensures that a given Kt0 , unlike the classic in-

506

dex Kt , will represent meteorologicaly similar conditions regardless of the sun’s

507

position.. Note that error statistics reported in Table 5 were calculated from

508

validation sites accounting for at least three residential PV installations (i.e.

509

the Middelkerke/AWS 6407, Melle/AWS 6434 and Stabroek/AWS 6438 mea-

510

surement stations were excluded). Clearly, the relative accuracy of the GHI

511

estimates varies noticeably as the sky conditions moves from overcast to clear

512

sky situations irrespective of the calculation method. With a reported RMSE

513

value of roughly 116%, the 1− PV-M and X− PV-EM approaches fail to pro-

514

duce reliable estimations in overcast conditions (i.e. 0.0 ≤ Kt0 < 0.2) where

515

the global radiation is mainly composed of diffuse radiation. With a reported

516

RMSE in the order of 70% the satellite-based retrieval method also exhibits

517

a poor performance in overcast conditions but shows a rapid performance im-

518

provement as the sky becomes clearer and presents a minimum RMSE value of

AC C

EP

TE D

M AN U

SC

502

lar irradiance, Gh , normalized by the corresponding extra-atmospheric irradi −1 −1.253 ance, E = Io cos θz ). M = sin γ + 0.15(γ + 3.885) is the optical

501

519

11.3% in partly clear conditions (i.e. 0.6 ≤ Kt0 < 0.8). Similarly, the accuracy

520

of GHI estimations from measured AC PV output power increases as the sky

521

conditions becomes clear but the magnitude of the errors is still at least twice

522

the one found for the satellite retrieval method irrespective of the computation

523

approaches. As an example, RMSE values of roughly 30% are reported for the

524

1− PV-M and X− PV-EM approaches and of about 17% for the satellite-based

525

method in the clear sky bin (i.e. 0.8 ≤ Kt0 ≤ 1.0) where the direct component 21

ACCEPTED MANUSCRIPT

largely domines. The poor performance of the PV systems power output method

527

in overcast and cloudy conditions (i.e. 0.0 ≤ Kt0 < 0.4) is mainly due to the

528

reduced accuracy in the derivation of the in-plane hourly global solar irradiance

529

from the hourly PV system energy production as computed from Eq. 17 when

530

the global radiation is dominated by diffuse radiation (i.e. at low PV module’s

531

yield). This is supported by the recent evaluation of the tilt to horizontal global

532

solar irradiance conversion by [29] indicating that the relative accuracy of the

533

conversion does not vary noticeably as the sky condition moved from overcast

534

to clear sky situations irrespective of the considered transposition model (e.g. a

535

RMSE variation of about 3.5% was reported for the SKA transposition model

536

we use in the current study).

537

M AN U

SC

RI PT

526

Finally, both the PV systems power output method and the MSG based

539

retrieval method present a positive bias when the global radiation is dominated

540

by diffuse radiation (with MBE values of e.g. 43.1%, 37.3% and 45.4% reported

541

in the overcast bin for the 1− PV-M approach, the X− PV-EM approach and the

542

MSG based method, respectively) and a positive bias under clear-sky situation

543

(with e.g. MBE values of respectively -18.6%, -13.% and -9.7% reported in the

544

clear-sky bin). In overall the MSG based retrieval method present a positive

545

bias in the order of 3.5% while the PV systems power output method shows

546

a negative bias for both the 1− PV-M and X− PV-EM approaches (i.e. overall

547

MBE of -3.4% and -1.%, respectively).

548

4.2. Daily validation

AC C

EP

TE D

538

549

Computation of the statistical errors indexes on a daily basis is not as

550

straightforward as for an hourly basis because as already mentioned both the

551

1− PV-M and X− PV-EM approaches fail to produce valid GHI estimates at low

552

solar elevation conditions. In Table 6 RMSE and MBE indexes have been com-

553

puted by assuming no incoming global horizontal solar irradiance in the com-

554

putation of the daily global horizontal solar irradiation for data points where

555

no valid hourly GHI estimates were obtained. Because unsuccessful GHI esti22

ACCEPTED MANUSCRIPT

mations from AC power output can be as large as 39.7% when only one PV

557

installation is considered (see Table 3) the daily performance of the PV method

558

can be very limited in some places in comparison to the satellite method. As

559

an example, the RMSE of 1321.2 W.h.m−2 or 42.9% and the MBE of -984.8

560

W.h.m−2 or -31.9% reported at the Middelkerke validation site (AWS 6407) in

561

Table 6 for the PV system method reduce to a RMSE of 296.4 W.h.m−2 or 9.6%

562

and a MBE of 43 W.h.m−2 or 1.4% for the MSG-based retrieval method at this

563

validation site.

SC

RI PT

556

564

Globally Table 6 indicates that the daily computation exhibits a better per-

566

formance than hourly estimation irrespectively of the considered retrieval meth-

567

ods. An overall RMSE of 11% and a slight positive bias is reported for the satel-

568

lite based method. 1− PV-M and X− PV-EM behave quite similarly in terms of

569

RMSE (i.e. overall RMSE of 12%). Clearly the RMSE magnitude difference

570

between the MSG-based and the PV systems power output methods is drasti-

571

cally reduced when considering daily global solar irradiation quantities rather

572

than hourly GHI values. However the RMSE spatial variation (i.e. from one

573

validation site to another) is larger in the PV systems based method and, while

574

a systematic positive bias is reported for the satellite retrieval method, the sign

575

of the bias can vary from one validation site to another in the PV systems based

576

method and even between the 1− PV-M and X− PV-EM approaches on a given

577

site.

TE D

EP

AC C

578

M AN U

565

579

4.3. Spatial validation

580

Finally, to assess the spatial distribution of the solar surface irradiation

581

computed from hourly PV systems power outputs, the 1470 residential PV in-

582

stallations considered in our study (see Figure 1) were spatially aggregated into

583

150 clusters (see black dots on panel C in Figure 3 for the clusters spatial dis-

584

tribution) using the k-means algorithm ([39]). K-means clustering partitions

585

a dataset into a small number of clusters by minimizing the distance between 23

ACCEPTED MANUSCRIPT

each data point and the center of the cluster it belongs to. A minimum of four

587

PV installations by clusters was imposed except for two of them located in the

588

vicinity of the Belgian coast because of the very low density of PV installations

589

found in this area.

RI PT

586

590

Figure 3 presents an example of daily solar surface irradiation over the Bel-

592

gian territory as computed from the interpolation of ground measurements us-

593

ing the ordinary kriging (OK) method (e.g. [40]), the MSG satellite derived

594

estimation and the PV systems method using X− PV-EM. Clearly, due to the

595

sparsity of the ground stations networks, interpolating ground data generates

596

only a coarse distribution of the solar surface irradiation: large-scale variations

597

of the solar irradiation (such as the south–east to north–west positive gradi-

598

ent) are identified but local fluctuations remain unseen (Figure 3, panel a).

599

Satellite-derived estimates, on the other hand, provide a global coverage and

600

are therefore able to account for clouds induced small-scale variability in sur-

601

face solar radiation (Figure 3, panel b). Regarding the PV systems method, the

602

daily solar irradiation estimated at the PV clusters level were then interpolated

603

by OK method to cover the entire Belgian territory (Figure 3, panel c). As we

604

see, the south–east to north–west positive gradient is well apparent as well as

605

some of the regional specificities. For instance, the Gaume region (area in the

606

south-east of Belgium) located on the south side of the Ardenne (hilly mass) and

607

that enjoys longer sunshine time appears clearly on the mapping. In general,

608

the PV systems method provides small-scale patterns partly supported by the

609

MSG derived mapping. Some others appear as the signature of an erroneous

AC C

EP

TE D

M AN U

SC

591

610

estimation at the cluster level.

611

612

5. Conclusions and perspective

613

The accuracy of solar irradiation data derived from the energy production of

614

PV systems depends on a number of factors including the efficiency of the PV

24

ACCEPTED MANUSCRIPT

system, the weather conditions, the density of PV systems that can be used for

616

the tilt to horizontal conversion, other data sources that can be accessed to com-

617

plement the PV data. In particular, because the solar irradiation data obtained

618

by reverse engineering of the PV energy output data degrade as the informa-

619

tion about the orientation and tilt angles of the PV generator becomes more

620

inaccurate, reliable angular characterization of the PV generator is of prime

621

importance. Unfortunately, it was found that the provided angles can typically

622

bear inaccuracies up to 5-10o . To overcome such a limitation, the development

623

of a procedure able to assess the true azimuth and tilt angle of a PV genera-

624

tor from the sole knowledge of its hourly energy output data would be largely

625

beneficial for the method. Towards this objective, it is worth mentioning that

626

while the identification of a PV system’s location and orientation from perfor-

627

mance data has received limited attention in the past, [41] have introduced a

628

method for automatic detection of PV system configuration. And even more

629

recently, [42] have proposed a method to estimate the location and orientation

630

of distributed photovoltaic systems from their generation output data.

SC

M AN U

TE D

631

RI PT

615

Another limitation of the method is that there are certain sun positions for

633

which the method fails to produce a valid estimation. As an example, unsuc-

634

cessful tilt to horizontal conversions systematically occurs at low solar elevation

635

(i.e. γ ≤ 5o ) irrespective of the number of PV systems involved in the conversion

636

process. Increasing the number of PV installations involved in the computation

637

process allows reducing the conversion failure rate and the smoothing the es-

638

timation to some extent. As an example, validation results computed on an

AC C

EP

632

639

hourly basis provide a mean RMSE value of about 40% when considering a

640

group of neighboring PV installations in the estimation process while values

641

as large as 60% are reported when using a single installation. By comparison,

642

satellite-based global horizontal irradiance estimation exhibits a better perfor-

643

mance with an associated overall RMSE of about 20%. By contrast the overall

644

performance of both methods does not differ significantly on a daily basis. How-

645

ever, spatially the distribution of the solar surface irradiation computed from 25

ACCEPTED MANUSCRIPT

PV systems outputs presents small-scale pattern artifacts signature of erroneous

647

estimations.

RI PT

646

648

Beside a more accurate characterization of the orientation and tilt angles

650

of the PV generator, the conversion of PV system energy production to tilted

651

global solar irradiance could certainly be improved by taking into account the

652

wind influence on the operating temperature of the solar cell as an example or by

653

performing clear-sky radiative transfer computations with a more sophisticated

654

model, etc. But because of the nature of the uncertainties that we have faced

655

during our study, it could also be valuable to investigate on how recent data min-

656

ing techniques (i.e. the extraction of hidden predictive information from large

657

databases [43]) can be incorporated in the transformation of electrical variables

658

into environmental variable such as solar irradiation. To this respect it is worth

659

pointing out that artificial neural networks and machine learning techniques

660

have already revealed to be powerful in the forecasting of solar irradiation data

661

(e.g. [44], [45], [46], [47] and [48]).

662

Acknowledgments

M AN U

TE D

663

SC

649

This study was supported by the Belgian Science Policy Office (BELSPO) through the Belgian Research Action through Interdisciplinary Networks (BRAIN-

665

be) pioneer research project BR/314/PI/SPIDER Solar Irradiation from the

666

energy production of residential PV systems .

EP

664

A. Decomposition models

668

A.1. ERB model (Erbs et al., 1982):

AC C 667

669

[22] developed a correlation between the hourly clearness index, Kt , and the

670

corresponding diffuse fraction, Kd based on 5 stations data. In each station,

671

hourly values of direct and global irradiances on a horizontal surface were regis-

672

tered. Diffuse irradiance was obtained as the difference of these quantities. The

26

ACCEPTED MANUSCRIPT

proposed correlation combines a linear regression for 0 < Kt ≤ 0.22, a fourth

674

degree polynom for 0.22 < Kt ≤ 0.8 and a constant value for Kt > 0.8

RI PT

673

   1 − 0.09 Kt    Kd = 0.9511 − 0.1604 Kt + 4.388 Kt2 − 16.638 Kt3 + 12.336 Kt4     0.165

Kt ≤ 0.22

0.22 < Kt ≤ 0.8 Kt > 0.8

675

A.2. OLS model (Skartveit and Olseth, 1987)

SC

(29)

Based on 15 years (i.e. 1965-1979) of hourly records of global and diffuse

677

horizontal irradiances with mean solar elevation larger than 10o performed in

678

Bergen (Norway: 60.4o N, 5.3o E) [49] proposed an analytical model expressing

679

the hourly diffuse fraction of global irradiance in terms of hourly solar elevation

680

and clearness index. In this model, the direct normal irradiance, Bn , is derived

681

from the global horizontal irradiance, Gh , and the solar elevation angle, γ, using

682

the following equation:

M AN U

676

Gh (1 − Ψ) sin γ

(30)

TE D

Bn = 683

where Ψ is a function of the clearness index, Kt . The model was validated with

684

data collected in 12 stations worldwide. The function Ψ reads as:

685

EP

   1 for Kt < c1    h i Ψ = 1 − (1 − d1 ) d2 c1/2 + (1 − d2 ) c23 for c1 ≤ Kt ≤ 1.09c2 3     1 − 1.09c 1−Υ for Kt > 1.09c2 2 Kt where:

687

c1 = 0.2

AC C 686

688

689

c2 = 0.87 − 0.56 −0.06 γ 
 c3 = 0.5 1 + sin π( dc43 − 0.5)

690

c4 = Kt − c1

691

d1 = 0.15 + 0.43 e−0.06 γ

692

d2 = 0.27

27

(31)

ACCEPTED MANUSCRIPT

694

695

d3 = c2 − c1   1/2 Υ = 1 − (1 − d1 ) d2 cc3 + (1 − d2 ) cc23 
 4 cc3 = 0.5 1 + sin π( cc d3 − 0.5)

696

cc4 = 1.09 c2 − c1

697

B. Transposition models

RI PT

693

Each model develops the diffuse transposition factor (i.e., the ratio of diffuse

699

radiation on a tilted surface to that of a horizontal), Rd , according to specific

700

assumptions.

701

B.1. HAY model (Hay, 1979):

M AN U

SC

698

702

In the HAY model ([23]), diffuse radiation from the sky is composed of an

703

isotropic component and a circumsolar one. Horizon brightening is not taken

704

into account. An anisotropy index, FHay , is used to quantify a portion of the

705

diffuse radiation treated as circumsolar with the remaining portion of diffuse

706

radiation assumed to be isotropic, i.e.,

TE D

RdHAY = FHAY rb + (1 − FHAY )

1 + cos β 2

(32)

where FHAY = Bh /(Io cos θz ) is the Hay’s sky-clarity factor and rb = cos θi / cos θz

708

the beam radiation conversion factor. θz is the solar zenith angle, θi is the in-

709

cidence angle of the beam radiation on the tilted surface, Bh is the direct hori-

710

zontal solar irradiance and Io is the extraterrestrial normal incident irradiance.

711

B.2. SKA model (Skartveit and Olseth, 1986):

EP

707

Solar radiation measurements indicate that a significant part of sky diffuse

713

radiation under overcast sky conditions comes from the sky region around the

714

zenith. This effect vanishes when cloud cover disappears. [50] modified the HAY

715

model (Eq. 32) in order to account for this effect,

AC C

712

RdSKA = FHAY rb + Z cos β + (1 − FHAY − Z)

1 + cos β 2

(33)

716

where Z = max(0, 0.3 − 2FHAY ) is the Skartveit-Olseth’s correction factor. If

717

FHAY ≥ 0.15, then Z = 0 and the model reduces to the HAY model. 28

ACCEPTED MANUSCRIPT

718

C. Description of the RMIs MAGIC/Heliosat-2 algorithm Information on GHI over Belgium is retrieved from MSG data by an adapta-

720

tion of the well-known Heliosat-2 method ([51]) initially developed for Meteosat

721

First generation satellites. The principle of the method is that a difference

722

in global radiation perceived by the sensor aboard a satellite is only due to a

723

change in the apparent albedo, which is itself due to an increase of the radiation

724

emitted by the atmosphere towards the sensor (i.e. [52], [53]). A key parameter

725

is the cloud index (also denoted as effective cloud albedo), n, determined by the

726

magnitude of change between what is observed by the sensor and what should

727

be observed under a very clear sky. To evaluate the all-sky GHI, a clear-sky

728

model is coupled with the retrieved cloud index which acts as a proxy for cloud

729

transmittance. Inputs to the Heliosat-2 method are not the visible satellite im-

730

ages in digital counts as in the original version of the method ([52]) but images

731

of radiances/reflectances:

M AN U

SC

RI PT

719

ρt (i, j) − ρtcs (i, j) ρtmax (i, j) − ρtcs (i, j)

(34)

TE D

nt (i, j) =

where nt (i, j) is the cloud index at time t for the satellite image pixel (i, j);

733

ρt (i, j) is the reflectance or apparent albedo observed by the sensor at time t;

734

ρtmax (i, j) is the apparent albedo of the brightest cloud at time t; ρtcs (i, j) is the

735

apparent ground albedo under clear-sky condition at time t.

736

EP

732

737

The MAGIC/Heliosat-2 method implemented at RMI ([54]) is applied to

738

the visible narrow-bands of the Spinning Enhanced Visible and Infrared Imager (SEVIRI) on board of the MSG platform. ρtcs (i, j) is derived from the fifth

740

percentile of the ρt (i, j) distribution related to the last 60 days ([55]) while

741

ρtmax (i, j) is estimated from theoretical radiative transfer model calculation (see

742

[54]). GHI estimates are computed for each pixel and MSG time slot by the

743

combination of the satellite cloud index, n, and the GHI in clear sky condition

744

calculated by the Mesoscale Atmospheric Global Irradiance Code (MAGIC)

745

clear-sky model ([56]). Daily global solar radiations are computed by trapezoidal

AC C 739

29

ACCEPTED MANUSCRIPT

integration over the diurnal cycle of the retrieved GHI t (i, j). The retrieval

747

process runs over a spatial domain ranging from 48.0o N to 54.0o N and from

748

2.0o E to 7.5o E within the MSG field–of–view. In this domain, the SEVIRI

749

spatial sampling distance degrades to about 6 km in the north-south direction

750

and 3.3 km in the east-west direction.

751

References

SC

RI PT

746

[1] Perez R, Seals R, Stewart R, Zelenka A, Estrada-Cajigal V. Using satellite-

753

derived insolation data for the site/time specific simulation of solar energy

754

systems. Solar Energy 1994;53:491–5.

M AN U

752

755

[2] Perez R, Seals R, Zelenka A. Comparing satellite remote sensing and ground

756

network measurements for the production of site/time specific irradiance

757

data. Solar Energy 1997;60:89–96.

[3] Hammer A, Heinemann D, Hoyer C, Kuhlemann R, Lorentz E, Mueller R,

759

et al. Solar energy assessment using remote sensing technologies. Remote

760

Sensing of Environment 2003;86:423–32.

TE D

758

[4] Renn´e D, George R, Wilcox S, Stoffel T, Myers D, Heimiller D. Solar re-

762

source assessment. technical report nrel/tp-581-42301. National Renewable

763

Energy Laboratory 2008;1617 Cole Boulevard, Golden, Colorado:80401–

764

3393.

EP

761

[5] Killinger S, Braam F, Mueller B, Wille-Haussmann B, McKenna R. Pro-

766

jection of power generation between differently-oriented pv systems. Solar

767

Energy 2016;136:153–65.

AC C

765

768

[6] Elsinga B, van Sark W, Ramakers L. Inverse photovoltaic yield model

769

for global horizontal irradiance reconstruction; vol. doi: 10.1002/ese3.162.

770

2017.

30

ACCEPTED MANUSCRIPT

[7] Leloux J, Narvarte L, Trebosc D. Review of the performance of residen-

772

tial pv systems in belgium. Renewable and Sustainable Energy Reviews

773

2012;16:178–84.

774

775

RI PT

771

[8] Journ´ee M, Bertrand C. Quality control of solar radiation data within the rmib solar measurements network. Solar Energy 2011;85:72–86.

[9] King D, Boyson W, Kratochvill J. Photovoltaic Array Performance Model.

777

Sandia National Laboratories (SANDIA REPORT SAND2004 3535): Al-

778

buquerque, New Mexico, 42 pp.; 2004.

SC

776

[10] Lorenzo E. Energy Collected and Delivered by PV modules. In Handbook

780

of Photovoltaic Science and Engineering, edited by Antonio Luque and

781

Steven Hegedus: John Wiley and Sons, Ltd, pp. 984–1042; 2011.

M AN U

779

782

[11] Decker B, Jahn U. Performance of 170 grid connected pv plants in northern

783

germany - analysis of yields and optimization potentials. Solar Energy

784

1997;59:127–33.

[12] Marion B, Adelstein J, Boyle K, Hayden H, Hammond B, Fletcher T, et al.

786

Performance parameters for grid-connected PV systems. In Conference

787

record of the Thirty-first IEEE Photovoltaic Specialists Conference: pp.

788

1601–1606; 2005.

TE D

785

[13] Mart´ınez-Moreno F, Lorenzo E, Narvarte L, Moreton R, Tyutyundzhiev

790

N, De La Parra I. Technical specifications and quality control procedures

791

for reducing the uncertainty in PV installations: results of the FP7 project

AC C

EP

789

792

PVCROPS. In Proceedings of the 31st European Photovoltaic Solar Energy

793

Conference and Exhibition (EU PVSEC): pp. 2217–2221; 2015.

794

795

796

797

[14] Evans D. Simplified method for predicting photovoltaic array output. Solar Energy 1981;27:555–60.

[15] Osterwald C. Translation of device performance measurements to reference conditions. Solar Cells 1986;18:269–79.

31

ACCEPTED MANUSCRIPT

[16] Randall J, Jacot J. Is am1.5 applicable in practice? modelling eight photo-

799

voltaic materials with respect to light intensity and two spectra. Renewable

800

Energy 2003;28:1851–1864.

RI PT

798

[17] Williams S, Betts T, Helf T, Gottschalg R, Beyer H, Infield D. Long-term

802

module performance based on realistic reporting conditions with consider-

803

ation to spectral effects. In 3rd World Conference on Photovoltaic Energy

804

Conversion: Osaka, Japan, pp. 1908–1911; 2003.

806

[18] Mart`ın N, Ruiz J. Annual angular reflection losses in pv modules. Prog Photovolt: Res Appli 2005;13:75–84.

M AN U

805

SC

801

807

[19] Taylor J, Leloux J, Everard A, Briggs J, Buckley A. Monitoring thou-

808

sands of distributed systems in the UK: Energy production and perfor-

809

mance. In 11th Photovoltaic Science Application and Technology Confer-

810

ence: PVSAT-11: University of Leeds: Leeds, UK, Wednesday 15th-Friday

811

17th April 2015: Conference Proceedings C97; 2015.

[20] Reno M, Hansen C, Stein J. Global horizontal Irradiance Clear Sky Mod-

813

els: Implementation and Analysis. Sandia National Laboratories (SANDIA

814

REPORT SAN2012-2389): Albuquerque, New Mexico, 68 pp.; 2012.

816

[21] Ineichen P, Perez R. A new airmass independent formulation for the linke turbidity coefficient. Solar Energy 2002;73:151–7.

EP

815

TE D

812

[22] Erbs D, Klein S, Duffie J. Estimation of the diffuse radiation fraction

818

for hourly, daily and monthly-averaged global radiation. Solar Energy

AC C

817

819

1982;28:293–302.

820

[23] Hay J. Study of shortwave radiation on non-horizontal surfaces. Rep. No.

821

79-12, Atmospheric Environment Service: Downsview, Ontario; 1979.

822

[24] Yang D. Solar radiation on inclined surfaces: corrections and benchmarks.

823

Solar Energy 2016;136:288–302.

32

ACCEPTED MANUSCRIPT

825

[25] Faiman D, Zemel A, Zangvil A. A method for monitoring insolation in remote regions. Solar Energy 1987;38:327–33.

RI PT

824

[26] Yang D, Dong Z, Nobre A, Khoo Y, Jirutitijaroen P, Walsh W. Evaluation

827

of transposition and decomposition models for converting global solar irra-

828

diance form tilted surface to horizontal in tropical regions. Solar Energy

829

2013;97:369–87.

SC

826

[27] Yang D, Ye Z, Nobre A, Du H, Walsh W, Lim L, et al. Bidirectional irradi-

831

ance transposition based on the perez model. Solar Energy 2014;110:768–

832

80.

M AN U

830

833

[28] Marion B. A model for deriving the direct normal and diffuse horizontal

834

irradiance from the global tilted irradiance. Solar Energy 2015;122:1037–46.

835

[29] Housmans C, Ipe A, Bertrand C. Tilt to horizontal global solar irradiance

836

conversion: An evaluation at high tilt angles and different orientations.

837

Renewable Energy 2017;113:1529–38.

[30] Gueymard CA. Direct and indirect uncertainties in the prediction of tilted

839

irradiance for solar engineering applications. Solar Energy 2009;83:432–44.

840

[31] Olmo F, Vida J, Foyo I, Castro-Diez Y, Alados-Arboledas L. Prediction

841

of global irradiance on inclined surfaces from horizontal global irradiance.

842

Energy 1999;24:689–704.

EP

TE D

838

[32] Grosan C, Abraham A. A new approach for solving nonlinear equations

844

systems. IEEE Trans Syst Man Cybern Part A: Syst Humans 2008;38:698–

AC C

843

845

846

714.

[33] Powell M. An efficient method for finding the minimum of a function

847

of several variables without calculating derivatives.

848

1964;7(2):155–62.

849

850

Computer Journal

[34] Acton F. Numerical Methods That Work, 1990, corrected edition. Mathematical Association of America, pp. 464-467: Washington; 1970. 33

ACCEPTED MANUSCRIPT

[35] Press W, Teukolsky S, Vetterling W, Flannery B. Numerical recipes in C.

852

The art of Scientific Computing. Second Edition (Reprinted with correc-

853

tions, 1997). Cambridge University Press: Cambrige, UK, 994 pp.; 1992.

RI PT

851

[36] Schmetz J, Pili P, Tjemkes S, Just D, Kerkmann J, Rota S, et al. An

855

introduction to meteosat second generation (msg). Bulletin of the American

856

Meteorological Society 2002;83:977–92.

SC

854

[37] Perez R, Ineichen P, Seals R, Zelenka A. Making full use of the clear-

858

ness index for parameterizing hourly insolation conditions. Solar Energy

859

1990;45:111–4.

860

[38] Kasten F.

M AN U

857

A simple parameterization of the pyrheliometric formula

861

for determining the linke turbidity factor. Meteorologische Rundschau

862

1980;33:124–7.

864

865

866

[39] Hair J, Black W, Babin B, Anderson R. Multivariate Data Analysis. Seventh ed. Prentice Hall; 2009.

[40] Wackernagel H. Multivariate geostatistics: an introduction with applica-

TE D

863

tions. Berlin:Springer-Verlag; 1995. [41] Williams M, Kerrigan S, Thornton A, Energy L. Automatic Detection

868

of PV System Configuration. Presented at the World Renewable Energy

869

Forum: 5 pp.; 2012.

EP

867

[42] Haghdadi N, Copper J, Bruce A, MacGill I. A method to estimate the

871

location and orientation of distributed photovoltaic systems from their gen-

AC C

870

872

eration output data. Renewable Energy 2017;108:390–400.

873

[43] Han J, Kamber M. Data Mining: Concepts and Techniques, 2nd ed. March

874

2006. ISBN 1-55860-901-6.: The Morgan Kaufmann Publishers; 2006.

875

[44] Sfetsos A, Coonick A. Univariate and multivariate forecasting of hourly

876

solar radiation with artificial intelligence techniques.

877

2000;68:169–78. 34

Solar Energy

ACCEPTED MANUSCRIPT

[45] Behrang M, Assareh E, Noghrehabadi A. The potential of different artificial

879

neural network (ann) techniques in daily global solar radiation modeling

880

based on meteorological data. Solar Energy 2010;84:1468–80.

RI PT

878

[46] Mellit A, Massi Pavan A. A 24-h forecast of solar irradiance using artificial

882

neural network: Application for performance prediction of a grid-connected

883

pv plant at trieste, italy. Solar Energy 2010;84:807–21.

884

885

SC

881

[47] Almeida M, Perpi¨ n´ an O, Narvarte L. Pv power forecast using a nonparametric pv model. Solar Energy 2015;115:354–68.

[48] Voyant C, Notton G, Kalogirou S, Nivet ML, Paoli C, Motte F, et al. Ma-

887

chine learning methods for solar radiation forecasting: A review. Renewable

888

Energy 2017;105:569–82.

890

891

892

893

894

[49] Skartveit A, Olseth J. A model for the diffuse fraction of hourly global radiation. Solar Energy 1987;38:271–4.

[50] Skartveit A, Olseth J. Modelling slope irradiance at high latitudes. Solar Energy 1986;36:333–44.

TE D

889

M AN U

886

[51] Rigollier C, Lef`evre M, Wald L. The method heliosat-2 for deriving shortwave solar radiation from satellite images. Solar Energy 2004;77:159–69. [52] Cano D, Monget J, Albuisson M, Guillar H Regas N, Wald L. A method for

896

the determination of the global solar radiation from meteorological satellite

897

data. Solar Energy 1986;56:207–12.

EP

895

[53] Raschke E, Gratzki A, Rieland M. Estimates of global radiation at the

899

ground from reduced data sets of international satellite cloud climatology

900

project. J Climate 1987;7:205–13.

AC C 898

901

[54] Demain C, Journ´ee M, Bertrand C. Sensitivity of the rmi’s magic/heliosat-

902

2 method to relevant input data.

903

2013;10:7–13.

35

Advances in Sciences and Research

ACCEPTED MANUSCRIPT

[55] Ipe A, Clerbaux N, Bertrand C, Dewitte S, Gonzalez L. Pixel-scale com-

905

posite top-of-the-atmosphere clear-sky reflectances for meteosat-7 visible

906

data. J Geophys Res 2003;108:148–227.

RI PT

904

[56] M¨ uller R, Matsoukas C, Gratzki A, Behr H, Hollmann R. The cm-saf oper-

908

ational scheme for the satellite based retrieval of solar surface irradiance-a

909

lut based eigenvector hybrid approach. Remote Sensing of the Environment

910

2009;113:1012–24.

AC C

EP

TE D

M AN U

SC

907

36

ACCEPTED MANUSCRIPT

1

913 914 915

2

916 917 918 919 920 921 922 923

AC C

EP

TE D

M AN U

924

3

Location of the selected 1470 residential PV systems within the Belgian territory. Red dots indicate the location of the RMI’s ground radiometric stations considered in this study . . . . . . . Definition of angles used as coordinates for an element of sky radiation to an inclined plane of tilt β . . . . . . . . . . . . . . . Comparison between the daily spatial distribution of surface solar global irradiation (in W.h.m−2 ) over Belgium as computed (a) by the interpolating ground measurements, (b) by the MSG satellite retrieval method and, (c) by the PV systems power output method. Illustrations are for the 20 June, 2014. Black dots in panel (a) indicate the location of the RMI’s in-situ measurements sites. Black dots in panel (c) indicate the location of the 150 clusters of PV systems . . . . . . . . . . . . . . . . . . . . . .

RI PT

912

List of Figures

SC

911

37

38

39

40

ACCEPTED MANUSCRIPT

2˚30'

3˚00'

3˚30'

4˚00'

4˚30'

5˚00'

5˚30'

6˚00'

6˚30' 51˚30'

RI PT

51˚30'

51˚00'

51˚00'

50˚30'

SC

50˚30'

N

km 49˚30' 2˚30'

0

50 3˚00'

3˚30'

M AN U

50˚00'

4˚00'

4˚30'

5˚00'

5˚30'

6˚00'

50˚00'

49˚30' 6˚30'

AC C

EP

TE D

Figure 1: Location of the selected 1470 residential PV systems within the Belgian territory. Red dots indicate the location of the RMI’s ground radiometric stations considered in this study

38

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

AC C

EP

TE D

Figure 2: Definition of angles used as coordinates for an element of sky radiation to an inclined plane of tilt β

39

ACCEPTED MANUSCRIPT

(a) AWS 20/06/2014 W.h.m−2

●

RI PT

7500

●

●

7000

●

6500

● ●

●

6000

●

5500

● ●

SC

5000 4500

●

●

4000

M AN U

3500 3000

●

(b) MSG

TE D

20/06/2014

2500

W.h.m−2 7500 7000 6500 6000 5500 5000 4500 4000

EP

3500 3000 2500

AC C

(c) PV systems

20/06/2014

W.h.m−2 7500 ●

●

●

7000

●

●

●

●

●

6500

● ● ●

●

●

●

6000

● ●

●

● ● ●

●

●

●

● ● ●

● ●

●

●

●

●

●

●

●

● ●

●

●

●

●

●

● ● ●

●

● ●

●

●

●

●

●

● ●

● ●

● ● ●● ● ● ● ● ● ●● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ● ●

●

●

●

●

● ● ● ●

●

●

●

●

●

● ●

● ●

●●

● ●

●

●

5000

●

●

● ● ●

4500

●

●

5500

●

●

● ●

●

●

●

4000

●

● ●

● ●

40

●

●

3500

● ●

●

●

● ●

● ●

● ●

●

●

3000

● ● ● ● ● ● ●

2500

Figure 3: Comparison between the daily spatial distribution of surface solar global irradiation (in W.h.m−2 ) over Belgium as computed (a) by the interpolating ground measurements, (b) by the MSG satellite retrieval method and, (c) by the PV systems power output method. Illustrations are for the 20 June, 2014. Black dots in panel (a) indicate the location of the RMI’s in-situ measurements sites. Black dots in panel (c) indicate the location of the 150

ACCEPTED MANUSCRIPT

Table 1: Distribution of the number of PV systems as a function of the orientation and tilt angle. <–East -80 -70

0 4 4 13 67 4 0 0 0 0

0 0 2 4 8 2 0 0 0 0

0 0 1 14 21 2 0 0 0 0

-60

-50

-40

-30

0 0 2 4 13 2 0 0 0 0

0 0 0 0 0 0 0 0 0 0

1 6 19 50 188 12 0 0 0 0

0 0 10 24 38 3 0 0 0 0

Deviation from South (o ) -20 -10 0 10 20 0 2 12 24 36 6 0 0 0 0

0 2 10 18 30 4 0 0 0 0

3 7 20 71 223 16 0 0 0 0

0 0 8 9 16 2 0 0 0 0

0 2 5 21 32 3 1 0 0 0

AC C

EP

TE D

M AN U

SC

0 10 20 30 40 50 60 70 80 90

-90

41

30

40

50

60

West –> 70 80

0 0 9 14 23 1 0 0 0 0

0 0 0 0 0 0 0 0 0 0

0 2 18 47 147 17 0 0 0 0

0 0 2 7 24 5 0 0 0 0

0 0 4 6 12 6 0 0 0 0

RI PT

Tilt angle (o )

0 0 3 6 15 1 0 0 0 0

90 0 0 0 0 0 0 0 0 0 0

ACCEPTED MANUSCRIPT

Table 2: Location of the ground stations involved in the RMI solar radiation monitoring network in which the global horizontal radiation is recorded. Lat. (o N) 51.198 50.905 50.976 51.326 51.076 50.798 50.096 50.583 51.222 50.194 50.040 50.916 49.621 50.512

Lon. (o E) 2.869 3.123 3.825 4.365 4.526 4.359 4.596 4.691 5.028 5.257 5.405 5.451 5.589 6.075

Alt. (m) 3.0 25.0 15.0 5.0 10.0 101.0 233.0 157.0 21.0 296.0 557.0 39.0 324.0 673.0

RI PT

Station Name MIDDELKERKE BEITEM MELLE STABROEK SINT-KATELIJNE-WAVER UCCLE DOURBES ERNAGE RETIE HUMAIN SAINT-HUBERT DIEPENBEEK BUZENOL MONT RIGI

AC C

EP

TE D

M AN U

SC

Station Code 6407 6414 6434 6438 6439 6447 6455 6459 6464 6472 6476 6477 6484 6494

42

ACCEPTED MANUSCRIPT

Hourly values 4184 4529 4477 4666 / 4698 4346 4573 / 4583 4469 / /

0 12? 4? 37? 0 5 3? 10 0

PV Method 1− PV-M X− PV-EM 36.74% / 12.89% 17.40% 39.74% / 22.29% / / / 16.69% 18.09% 13.71% 15.55% 13.21% 21.21% / / 17.72% 23.04% 14.03% 22.04% / / / /

SC

PV systems number 1 5 1 1?

AC C

EP

TE D

M AN U

AWS code 6407 6414 6434 6438 6439 6447 6455 6459 6464 6472 6477 6484 6494

RI PT

Table 3: Unsuccessful conversion (in %) reported for the PV systems 1− PV-M and X− PV-EM approaches, respectively. Also provided is the total number of hourly data points available at each validation sites and the number of PV installations found in the vicinity of the measurement stations. ? indicates the PV installations located within a radius of 10 km surrounding the validation site.

43

6407 6414 6434 6438 6447 6455 6459 6472 6477

AWS code

5 3?

1 5 1 1? 12? 4? 37?

EP

0.54 8.68 11.61 13.86 7.67 15.63 1.01 11.26 13.32

MBE (0.16%) (3.11%) (3.88%) (5.01%) (2.84%) (5.87%) (0.36%) (4.18%) (4.80%)

Table 4: Comparison between global horizontal solar irradiance produced by the PV systems power output method for both the 1− PV-M and X− PV-EM approaches and retrieved from the MSG satellite images with the corresponding ground measurements. The RMSE and MBE error statistics (in W.m−2 and %) are calculated on a hourly basis over the full year 2014. ? indicates the PV installations located within a radius of 10 km surrounding the validation site.

RI PT

MSG RMSE (21.22%) (20.28%) (20.93%) (22.69%) (19.64%) (22.37%) (20.14%) (19.66%) (19.84%)

SC

69.55 56.68 62.58 62.80 53.13 59.52 57.00 52.95 55.05

M AN U

X− PV-EM approach RMSE MBE / / / / 113.81 (40.73%) -9.00 (-3.22%) / / / / / / / / 107.76 (39.83%) -9.03 (-3.34%) 120.87 (45.42%) 14.63 (5.50%) 129.64 (45.81%) -15.84 (-5.60%) 129.87 (48.22%) 7.86 (2.92%) 129.15 (46.55%) -3.56 (-1.28%)

TE D

PV system 1− PV-M approach RMSE MBE 189.34 (57.76%) -53.30 (-16.26%) 116.83 (41.81%) -24.06 (-8.61%) 116.94 (39.11%) -10.48 (-3.50%) 128.05 (46.26%) -32.51 (-11.75%) 106.14 (39.23%) 1.85 (0.68%) 114.88 (43.17%) 20.86 (7.84%) 112.23 (39.66%) -28.03 (-9.90%) 120.55 (44.76%) -9.08 (-3.37%) 110.13 (39.70%) -17.95 (-6.47%)

AC C

Nbr

ACCEPTED MANUSCRIPT

44

ACCEPTED MANUSCRIPT

73.83 106.96 117.59 133.60 128.89

(117.00%) (64.79%) (38.78%) (28.56%) (30.09%)

Method X− PV-EM 72.70 (115.21%) 112.55 (68.17%) 133.27 (43.95%) 142.86 (30.54%) 134.70 (31.45%)

113.49

(41.37%)

121.87

1− PV-M

(44.43%)

MSG

Data Number

46.26 54.82 63.48 53.04 72.45

(73.31%) (33.20%) (20.94%) (11.34%) (16.91%)

4164 5494 4981 6214 863

(19.18%) (25.30%) (22.94%) (28.61%) (3.97%)

55.75

(20.32%)

21716

(100.0%)

AC C

EP

TE D

M AN U

SC

Sky Condition 0.0 ≤ Kt0 < 0.2 0.2 ≤ Kt0 < 0.4 0.4 ≤ Kt0 < 0.6 0.6 ≤ Kt0 < 0.8 0.8 ≤ Kt0 ≤ 1.0 All Sky

RI PT

Table 5: Performance in term of RMSE of the PV systems power output method for both the 1− PV-M and X− PV-EM approaches and the MSG retrieval method in the global horizontal solar irradiance estimation as a function of the sky condition. Absolute (in W.m−2 ) and relative (in %) RMSE are computed from 21716 hourly 2014 data points. Validation sites accounting for less than 3 residential PV installations were discarded for the errors indices computation (i.e. AWS 6407, AWS 6434 and AWS 6438).

45

6407 6414 6434 6438 6447 6455 6459 6472 6477

AWS code

5 3?

EP

approach MBE / / -106.57 (-4.92%) / / / / -115.75 (-3.92%) 116.07 (6.56%) -64.94 (-3.02%) 43.76 (2.74%) 100.60 (3.98%)

MSG MBE (1.40%) (5.37%) (4.60%) (5.31%) (3.19%) (8.12%) (1.90%) (7.04%) (6.37%)

Table 6: Comparison between daily cumulated surface solar irradiation produced by the PV systems power output method (for both the 1− PV-M and X− PV-EM approaches) and retrieved from the MSG satellite images with the corresponding daily ground measurements. The RMSE and MBE error statistics (in W.h.m−2 and %) are calculated on a daily basis over the full year 2014.? indicates the PV installations located within a radius of 10 km surrounding the validation site.

RI PT

43.03 116.29 150.23 154.50 94.28 143.62 40.94 112.68 161.13

SC

RMSE 296.38 (9.61%) 268.53 (12.41%) 337.75 (10.35%) 322.61 (11.09%) 232.35 (7.87%) 260.73 (14.74%) 231.97 (10.79%) 185.28 (11.58%) 271.98 (10.76%)

M AN U

X− PV-EM RMSE / / 304.64 (14.08%) / / / / 242.89 (8.23%) 306.30 (17.32%) 180.90 (8.42%) 226.74 (14.18%) 259.35 (10.26%)

TE D

PV systems 1− PV-M approach RMSE MBE 1321.23 (42.85%) -984.80 (-31.94%) 308.40 (14.25%) -175.35 (-8.10%) 534.81 (16.39%) -257.06 (-7.88%) 640.95 (22.04%) -434.83 (-14.95%) 250.26 (8.48%) 5.98 (0.20%) 348.20 (19.69%) 168.84 (9.55%) 225.21 (10.48%) -137.49 (-6.40%) 151.70 (9.48%) -13.18 (-0.82%) 266.96 (10.56%) -29.97 (-1.19%)

AC C

1 5 1 1? 12? 4? 37?

Nbr

ACCEPTED MANUSCRIPT

46