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

     1  package eeslism
     2  
     3  import (
     4  	"fmt"
     5  	"math"
     6  	"strconv"
     7  	"strings"
     8  	"unicode"
     9  )
    10  
    11  /* VCFILE からの気象データ入力 */
    12  
    13  func wdflinit(Simc *SIMCONTL, Estl *ESTL, Tlist []TLIST) {
    14  	var wp WDPT
    15  	var s, ss, Err string
    16  	var dt float64
    17  	var id, N, i, m int
    18  
    19  	if s = Estl.Wdloc; s == "" {
    20  		return
    21  	}
    22  
    23  	Err = fmt.Sprintf(ERRFMT, "(wdflinit)")
    24  	Simc.Loc = NewLOCAT()
    25  	loc := Simc.Loc
    26  
    27  	m = -1
    28  
    29  	for _, field := range strings.Fields(s) {
    30  		N = len(field)
    31  		if st := strings.Index(field, "="); st >= 0 {
    32  			ss = field[:st]
    33  			dt, _ = strconv.ParseFloat(field[st+1:], 64)
    34  			switch ss {
    35  			case "Lat":
    36  				loc.Lat = dt
    37  			case "Lon":
    38  				loc.Lon = dt
    39  			case "Ls":
    40  				loc.Ls = dt
    41  			case "Tgrav":
    42  				loc.Tgrav = dt
    43  			case "DTgr":
    44  				loc.DTgr = dt
    45  			case "daymx":
    46  				loc.Daymxert = int(dt)
    47  			default:
    48  				id = 1
    49  			}
    50  		} else {
    51  			if field == "-" || m >= 0 {
    52  				switch field {
    53  				case "-Twsup":
    54  					m = 0
    55  				default:
    56  					loc.Twsup[m], _ = strconv.ParseFloat(field, 64)
    57  					m++
    58  					if m > 11 {
    59  						m = -1
    60  					}
    61  				}
    62  			} else {
    63  				loc.Name = field
    64  			}
    65  		}
    66  		s = s[N:]
    67  		for len(s) > 0 && unicode.IsSpace(rune(s[0])) {
    68  			s = s[1:]
    69  		}
    70  	}
    71  
    72  	if id != 0 {
    73  		fmt.Printf("%s %s\n", Err, ss)
    74  	}
    75  
    76  	wp.Ta = nil
    77  	wp.Xa = nil
    78  	wp.Rn = nil
    79  	wp.Ihor = nil
    80  	wp.Rh = nil
    81  	wp.Cc = nil
    82  	wp.Wv = nil
    83  	wp.Wdre = nil
    84  
    85  	for i = 0; i < Estl.Ndata; i++ {
    86  		t := &Tlist[i]
    87  		if t.Name == "Wd" {
    88  			s = t.Id
    89  			val := t.Fval
    90  			switch s {
    91  			case "T": // 温度
    92  				wp.Ta = val
    93  			case "x": // 絶対湿度
    94  				wp.Xa = val
    95  			case "Idn": // 法線面直達日射量
    96  				wp.Idn = val
    97  			case "Isky": // 水平面天空日射量
    98  				wp.Isky = val
    99  			case "Ihor": // 水平面全天日射量
   100  				wp.Ihor = val
   101  			case "CC": // 雲量
   102  				wp.Cc = val
   103  			case "Wdre": // 風向
   104  				wp.Wdre = val
   105  			case "Wv": // 風速
   106  				wp.Wv = val
   107  			case "RH": // 相対湿度
   108  				wp.Rh = val
   109  			case "RN": // 夜間放射量
   110  				wp.Rn = val
   111  			}
   112  		}
   113  	}
   114  
   115  	Simc.Wdpt = wp
   116  }
   117  
   118  func Wdflinput(wp *WDPT, Wd *WDAT) {
   119  	var Br float64
   120  
   121  	Wd.T = wp.Ta[0]
   122  	Wd.Idn = wp.Idn[0]
   123  	Wd.Isky = wp.Isky[0]
   124  
   125  	if wp.Ihor == nil {
   126  		Wd.Ihor = Wd.Idn*Wd.Sh + Wd.Isky
   127  	}
   128  
   129  	if wp.Xa != nil {
   130  		Wd.X = wp.Xa[0]
   131  	} else {
   132  		Wd.X = -999.0
   133  	}
   134  
   135  	if wp.Rh != nil {
   136  		Wd.RH = wp.Rh[0]
   137  	} else {
   138  		Wd.RH = -999.0
   139  	}
   140  
   141  	if wp.Cc != nil {
   142  		Wd.CC = wp.Cc[0]
   143  	} else {
   144  		Wd.CC = -999.0
   145  	}
   146  
   147  	if wp.Rn != nil {
   148  		Wd.RN = wp.Rn[0]
   149  	} else {
   150  		Wd.RN = -999.0
   151  	}
   152  
   153  	if wp.Wv != nil {
   154  		Wd.Wv = wp.Wv[0]
   155  	} else {
   156  		Wd.Wv = -999.0
   157  	}
   158  
   159  	if wp.Wdre != nil {
   160  		Wd.Wdre = wp.Wdre[0]
   161  	} else {
   162  		Wd.Wdre = 0.0
   163  	}
   164  
   165  	if Wd.X > 0.0 && Wd.RH < 0.0 {
   166  		Wd.RH = FNRhtx(Wd.T, Wd.X)
   167  	} else if Wd.X < 0.0 && Wd.RH > 0.0 {
   168  		Wd.X = FNXtr(Wd.T, Wd.RH)
   169  	}
   170  
   171  	if Wd.X > 0.0 {
   172  		Wd.H = FNH(Wd.T, Wd.X)
   173  	}
   174  
   175  	if Wd.X > 0.0 && Wd.CC > 0.0 || Wd.RN < 0.0 {
   176  		Br = 0.51 + 0.209*math.Sqrt(FNPwx(Wd.X))
   177  		Wd.RN = (1.0 - 0.62*Wd.CC/10.0) * (1.0 - Br) * Sgm * math.Pow(Wd.T+273.15, 4.0)
   178  		Wd.Rsky = ((1.0-0.62*Wd.CC/10.0)*Br + 0.62*Wd.CC/10.0) * Sgm * math.Pow(Wd.T+273.15, 4.0)
   179  	} else {
   180  		Wd.Rsky = Sgm*math.Pow(Wd.T+273.15, 4.0) - Wd.RN
   181  	}
   182  }