github.com/archlabjp/eeslism-go@v0.0.0-20231109122333-4bb7bfcdf292/eeslism/u_sun.go (about)

     1  /* ================================================================
     2  
     3   SUNLIB
     4  
     5    太陽位置および日射量計算用ライブラリ－
     6    （宇田川、パソコンによる空気調和計算法、プログラム4.1の C 言語版, ANSI C 版）
     7  
     8  ---------------------------------------------------------------- */
     9  
    10  package eeslism
    11  
    12  import (
    13  	"math"
    14  )
    15  
    16  func Sunint() {
    17  	var Rd float64 = math.Pi / 180.0
    18  	Slat = math.Sin(Lat * Rd)
    19  	Clat = math.Cos(Lat * Rd)
    20  	Tlat = math.Tan(Lat * Rd)
    21  	if UNIT == "SI" {
    22  		Isc = 1370.0
    23  	} else {
    24  		Isc = 1178.0
    25  	}
    26  }
    27  
    28  func FNDecl(N int) float64 {
    29  	return math.Asin(0.397949 * math.Sin(2.0*math.Pi*(float64(N)-81.0)/365.0))
    30  }
    31  
    32  func FNE(N int) float64 {
    33  	var B float64 = 2.0 * math.Pi * (float64(N) - 81.0) / 365.0
    34  	return 0.1645*math.Sin(2.0*B) - 0.1255*math.Cos(B) - 0.025*math.Sin(B)
    35  }
    36  
    37  func FNSro(N int) float64 {
    38  	return Isc * (1.0 + 0.033*math.Cos(2.0*math.Pi*float64(N)/365.0))
    39  }
    40  
    41  func FNTtas(Tt float64, E float64) float64 {
    42  	return Tt + E + (Lon-Ls)/15.0
    43  }
    44  
    45  func FNTt(Ttas float64, E float64) float64 {
    46  	return Ttas - E - (Lon-Ls)/15.0
    47  }
    48  
    49  func FNTtd(Decl float64) float64 {
    50  	var Tlat float64
    51  	var Cws, Ttd float64
    52  	Cws = -Tlat * math.Tan(Decl)
    53  	if 1.0 > Cws && Cws > -1.0 {
    54  		Ttd = 7.6394 * math.Acos(Cws)
    55  	} else {
    56  		if Cws >= 1.0 {
    57  			Ttd = 0.0
    58  		} else {
    59  			Ttd = 24.0
    60  		}
    61  	}
    62  	return Ttd
    63  }
    64  
    65  var __Solpos_Sdecl, __Solpos_Sld, __Solpos_Cld float64
    66  var __Solpos_Ttprev float64 = 25.0
    67  
    68  func Solpos(Ttas float64, Decl float64) (Sh float64, Sw float64, Ss float64, solh float64, solA float64) {
    69  	const PI float64 = math.Pi
    70  	var Ch, Ca, Sa, W float64
    71  
    72  	if Ttas < __Solpos_Ttprev {
    73  		__Solpos_Sdecl = math.Sin(Decl)
    74  		__Solpos_Sld = Slat * __Solpos_Sdecl
    75  		__Solpos_Cld = Clat * math.Cos(Decl)
    76  	}
    77  
    78  	W = (Ttas - 12.0) * 0.2618
    79  	Sh = __Solpos_Sld + __Solpos_Cld*math.Cos(W)
    80  	solh = math.Asin(Sh) / PI * 180.0
    81  
    82  	if Sh > 0.0 {
    83  		Ch = math.Sqrt(1.0 - Sh*Sh)
    84  		Ca = (Sh*Slat - __Solpos_Sdecl) / (Ch * Clat)
    85  		var fW0 float64
    86  		if W > 0.0 {
    87  			fW0 = 1.0
    88  		} else {
    89  			fW0 = -1.0
    90  		}
    91  		solA = fW0 * math.Acos(Ca) / PI * 180.0
    92  		Sa = (W / math.Abs(W)) * math.Sqrt(1.0-Ca*Ca)
    93  		Sw = Ch * Sa
    94  		Ss = Ch * Ca
    95  	} else {
    96  		Sh = 0.0
    97  		Sw = 0.0
    98  		Ss = 0.0
    99  		solh = 0.0
   100  		solA = 0.0
   101  	}
   102  
   103  	__Solpos_Ttprev = Ttas
   104  
   105  	return Sh, Sw, Ss, solh, solA
   106  }
   107  
   108  func Srdclr(Io float64, P float64, Sh float64, Idn *float64, Isky *float64) {
   109  	if Sh > 0.001 {
   110  		*Idn = Io * math.Pow(P, 1.0/Sh)
   111  		*Isky = Sh * (Io - *Idn) * (0.66 - 0.32*Sh) * (0.5 + (0.4-0.3*P)*Sh)
   112  	} else {
   113  		*Idn = 0.0
   114  		*Isky = 0.0
   115  	}
   116  }
   117  
   118  func Dnsky(Io float64, Ihol float64, Sh float64, Idn *float64, Isky *float64) {
   119  	if Sh > 0.001 {
   120  		Kt := Ihol / (Io * Sh)
   121  		if Kt >= 0.5163+(0.333+0.00803*Sh)*Sh {
   122  			*Idn = (-0.43 + 1.43*Kt) * Io
   123  		} else {
   124  			*Idn = (2.277 + (-1.258+0.2396*Sh)*Sh) * (Kt * Kt * Kt) * Io
   125  		}
   126  		*Isky = Ihol - *Idn*Sh
   127  	} else {
   128  		*Idn = 0.0
   129  		*Isky = Ihol
   130  	}
   131  }