import get from 'lodash-es/get.js'
import each from 'lodash-es/each.js'
import isNumber from 'lodash-es/isNumber.js'
import isearr from 'wsemi/src/isearr.mjs'
import isestr from 'wsemi/src/isestr.mjs'
import isnum from 'wsemi/src/isnum.mjs'
import isfun from 'wsemi/src/isfun.mjs'
import cdbl from 'wsemi/src/cdbl.mjs'
import _kriging from './interp2Kriging.mjs'
import aggregatePoints from './aggregatePoints.mjs'
import interp2Grid from './interp2Grid.mjs'
/**
* 不規則點資料先聚合再內插至規則網格
*
* Unit Test: {@link https://github.com/yuda-lyu/w-gis/blob/master/test/interp2Raster.test.mjs Github}
* @memberOf w-gis
* @param {Array} ops 輸入二維座標加觀測數據點陣列,為[{x:x1,y:y1,z:z1},{x:x2,y:y2,z:z2},...]點物件之陣列
* @param {Object} [opt={}] 輸入設定物件,預設{}
* @param {String} [opt.keyX='x'] 輸入點物件之x欄位字串,為座標,預設'x'
* @param {String} [opt.keyY='y'] 輸入點物件之y欄位字串,為座標,預設'y'
* @param {String} [opt.keyZ='z'] 輸入點物件之z欄位字串,為觀測值,預設'z'
* @param {String} [opt.modePick='min'] 輸入挑選方式字串,可選'min'、'max',預設'min'
* @param {Number} [opt.dx=0.000979087] 輸入網格x向間距數字,預設0.000979087,基於WGS84經緯度之數據點,代表約100m
* @param {Number} [opt.dy=0.000906247] 輸入網格y向間距數字,預設0.000906247,基於WGS84經緯度之數據點,代表約100m
* @param {Number} [opt.dxAgr=opt.dx*2] 輸入先聚合時網格x向間距數字,預設opt.dx*2
* @param {Number} [opt.dyAgr=opt.dy*2] 輸入先聚合時網格y向間距數字,預設opt.dy*2
* @param {Function} [opt.funValid=(x,y)=>{return true}] 輸入確認點座標(x,y)是否有效函數,回傳布林值,可使用Promise回傳,預設(x,y)=>{return true}
* @param {Function} [opt.funKriging=interp2Kriging] 輸入克利金處理函數,預設使用內建interp2Kriging
* @param {Function} [opt.funAdjust=(x,y,z)=>{return z}] 輸入內插後值調整函數,用於修正不合理值或做後處理,回傳布林值,可使用Promise回傳,預設(x,y,z)=>{return z}
* @param {Number} [opt.scale=1] 輸入當funKriging使用interp2Kriging時,正規化範圍數值,因處理多邊形時有數值容許誤差,故須通過縮放值域來減少問題,預設1是正規化0至1之間,使用scaleXY則是正規化為0至scaleXY之間,預設1
* @param {String} [opt.method='friedrich'] 輸入當funKriging使用interp2Kriging時,後端引擎字串,可選'friedrich'(Gaussian Process Regression,核函數)、'invall'(經典Ordinary Kriging,半變異圖擬合),預設'friedrich'
* @param {String} [opt.model='exponential'] 輸入當funKriging使用interp2Kriging時,擬合模式字串,method='invall'時可選'exponential'、'gaussian'、'spherical';method='friedrich'時'gaussian'對應SquaredExp、其餘(含spherical)對應Exponential,預設'exponential'
* @param {Number} [opt.sigma2=0] 輸入當funKriging使用interp2Kriging時,自動擬合參數sigma2數值,僅method='invall'時有效,預設0
* @param {Number} [opt.alpha=100] 輸入當funKriging使用interp2Kriging時,自動擬合參數alpha數值,僅method='invall'時有效,預設100
* @param {Number} [opt.ls=0.5] 輸入當funKriging使用interp2Kriging時,friedrich核函數lengthscale數值(對應正規化後值域),僅method='friedrich'時有效,預設0.5
* @param {Number} [opt.noise=0.01] 輸入當funKriging使用interp2Kriging時,friedrich觀測噪音標準差(對應正規化後值域),僅method='friedrich'時有效,預設0.01
* @returns {Object} 回傳規則網格物件,物件內grds為規則網格之二維陣列
* @example
*
* let ps
* let r
*
* ps = [
* { x: 121.48571681415929, y: 25.058614130434783, z: 95 },
* { x: 121.48500884955752, y: 25.062228260869563, z: 146 },
* { x: 121.47, y: 25.08125, z: 22 },
* { x: 121.53598230088495, y: 25.02475543478261, z: 223 },
* { x: 121.54582300884955, y: 25.074592391304346, z: 39 },
* { x: 121.52118584070796, y: 25.057472826086958, z: 191 },
* { x: 121.48083185840707, y: 25.02, z: 22 },
* { x: 121.50653097345132, y: 25.089429347826087, z: 249 },
* { x: 121.55, y: 25.09, z: 86 },
* { x: 121.52614159292035, y: 25.067364130434783, z: 125 }
* ]
* r = await interp2Raster(ps, { dx: 0.000979087 * 10, dy: 0.000906247 * 5 })
* console.log(JSON.stringify(r))
* // => {"xnum":9,"xmin":121.47,"xmax":121.55,"dx":0.00979087,"ynum":16,"ymin":25.02,"ymax":25.09,"dy":0.004531235,"zmin":21.97238320338241,"zmax":231.37651396132563,"grds":[[40.24080873716694,83.01239405163139,138.22453819053254,199.21833879316065,231.37651396132563,179.4092418651924,133.5966489628609,98.43131094726549,81.6157254369131],[28.612659600934037,78.64589776136685,134.92602500122075,190.68164433852962,213.40394727458326,174.83492333510154,129.08766696494985,89.8017417344036,67.72112466368993],[27.353273441923474,76.10407906727625,129.70239363820528,177.75803393621848,195.39203561609946,167.19063091217615,123.52852150481242,80.88069192638437,53.04342496772884],[35.72079384454732,75.77302739821397,124.04787860734024,165.59835724154416,181.6958063225275,160.28480895483338,118.89428840874089,76.92692555196564,43.935127229263585],[41.76867429377623,76.26387316623804,118.81264447813922,155.71687166433097,172.82503452816562,157.763181852382,117.24797520569288,82.26487092771593,54.69203567673017],[45.78653949445653,76.45547397137864,114.17040164627124,148.06086926717313,168.215091967896,163.3683127289209,125.90476566360043,94.41448179585412,70.05169881928228],[48.01594897112139,75.81223021661773,109.85625787409032,141.9070578353577,166.04010283449375,176.93358417480647,143.72203927754867,109.37409986509063,85.35057633465354],[48.65546052045457,73.98953210738432,105.2337168924808,136.35975170010016,163.6680042767658,186.22476224542035,157.94279936927273,124.3883692130921,99.99592081479383],[47.88516657866743,70.94060417873453,100.01423370698842,130.74693708572534,159.57133876509027,180.23089947582827,165.76541576119072,138.12446724788674,113.77324716229168],[45.87875479281789,66.82912494725846,94.28095114648264,124.78179507934202,154.19076671814713,174.5764477233918,170.72763187125065,150.61998513570094,126.6354125870738],[42.80196031075884,61.71254364579427,87.96979292064754,118.45484465132714,148.4684177959958,170.81268446853952,175.7048549274094,162.6215801592757,138.59900140619663],[38.8280984259269,55.586948934294156,81.04580917179617,111.86435292951454,142.84571515438563,168.38780143940264,181.64074496465543,174.9495861193142,149.53903286429568],[34.18939738163876,48.42874342247074,73.56680049388629,105.14028026991367,137.30016344590013,166.36896610133712,188.30452773986423,188.17294731217387,158.88604543814418],[29.283420811458402,40.210025965013415,65.78630667553657,98.46117276811188,131.6023672593224,163.56057128616015,193.89193971379814,202.11304791244845,165.3052258406751],[24.846404838076825,30.936935705360924,58.41581495379795,92.08453048361456,125.52559782167987,158.7990028603251,193.68886300636396,209.9813260869589,166.9018659615177],[21.97238320338241,22.075876720886487,52.89204643860076,86.29589538194874,118.98075936386894,151.54951416462308,184.075743199371,195.33462327956147,162.7980226006054]],"pts":[{"x":121.48571681415929,"y":25.058614130434783,"z":95},{"x":121.47,"y":25.08125,"z":22},{"x":121.53598230088495,"y":25.02475543478261,"z":223},{"x":121.54582300884955,"y":25.074592391304346,"z":39},{"x":121.52118584070796,"y":25.057472826086958,"z":191},{"x":121.48083185840707,"y":25.02,"z":22},{"x":121.50653097345132,"y":25.089429347826087,"z":249},{"x":121.55,"y":25.09,"z":86},{"x":121.52614159292035,"y":25.067364130434783,"z":125}]}
*
* ps = [
* { x: 121.48571681415929, y: 25.058614130434783, z: 95 },
* { x: 121.48500884955752, y: 25.062228260869563, z: 146 },
* { x: 121.47, y: 25.08125, z: 22 },
* { x: 121.53598230088495, y: 25.02475543478261, z: 223 },
* { x: 121.54582300884955, y: 25.074592391304346, z: 39 },
* { x: 121.52118584070796, y: 25.057472826086958, z: 191 },
* { x: 121.48083185840707, y: 25.02, z: 22 },
* { x: 121.50653097345132, y: 25.089429347826087, z: 249 },
* { x: 121.55, y: 25.09, z: 86 },
* { x: 121.52614159292035, y: 25.067364130434783, z: 125 }
* ]
* r = await interp2Raster(ps, { dx: 0.000979087 * 10, dy: 0.000906247 * 5, dxAgr: 0.000979087 * 15, dyAgr: 0.000906247 * 7.5 })
* console.log(JSON.stringify(r))
* // => {"xnum":9,"xmin":121.47,"xmax":121.55,"dx":0.00979087,"ynum":16,"ymin":25.02,"ymax":25.09,"dy":0.004531235,"zmin":22.109253252645797,"zmax":231.96077893080354,"grds":[[40.71659658523934,89.17668826319866,146.63335072669418,204.2701988704215,231.96077893080354,179.48157195765387,133.3487046345723,98.19979656600891,81.58221115212473],[28.890271410402164,88.33039035635299,147.44134632498577,198.7677236117231,216.21673647139988,175.61242305841537,129.04401761393652,89.63606532722903,67.69053894713387],[30.925335348278914,91.10830125096261,147.44059692517916,189.14084231100267,200.23199445044517,168.66016161963947,123.62905564209949,80.74465159848714,53.024008922962686],[45.34011643664557,97.7138752699657,147.86359597215778,179.77725106019355,187.93158645610606,162.20620917900385,119.00301962655178,76.76110080612706,43.90667151317203],[55.67253724493733,105.51822220598382,148.68027132615464,171.27469601591739,179.59734472101096,159.77044010196119,117.19010260566847,82.00578577656893,54.55952301452594],[61.316098206438156,110.02338063265654,146.669003186955,162.60612630606767,174.58289148793392,165.10004541507308,125.65933514656336,94.04398376451469,69.81196035551322],[62.17921631650031,101.53012771491981,130.03089022135083,152.91895841249595,171.22889523584595,178.08727742723303,143.3516594102452,108.91652532226541,85.03148171738144],[59.537840113738504,85.02829790162487,108.83472094581788,142.90131093411495,167.30790472727836,186.78962587366962,157.46911271859082,123.88427119316528,99.6301724090092],[55.32984680006676,75.22928886777967,99.98938320642368,133.7709306090002,161.72754167248584,180.53595440395745,165.277447091207,137.6164749683483,113.39300830696384],[50.69323186511775,68.65800379324453,93.55878843443456,125.71935665925217,155.18192307586364,174.6285564830647,170.26536170542587,150.14434369728062,126.269206596812],[45.86553597943461,62.49756374540951,87.0855265013598,118.33007437528377,148.66687693046322,170.64065806615537,175.27462979438465,162.20553468570674,138.26987709205997],[40.78693784256393,55.881749201388864,80.16290167186719,111.25811324787777,142.56732569148585,168.05475252285567,181.25043378916007,174.61428398429695,149.26341086007227],[35.464278345961716,48.49661904618487,72.7570825096106,104.35814448914948,136.77268473663332,165.94382743174845,187.96785334090742,187.93668661143275,158.6721052495077],[30.138119277026334,40.195220706639866,65.08507237219025,97.65815185220012,130.9733194967505,163.09923313157748,193.6176896965728,201.99019014876967,165.14993785248433],[25.44111089564515,30.926113328554862,57.83386031171706,91.33376283229556,124.88195272388997,158.33720416197545,193.4612406500904,209.94823508284605,166.78897000311312],[22.396555990499,22.109253252645797,52.41514593737879,85.62510783244069,118.3683319051399,151.1052849649574,183.85729392441013,195.27849783709144,162.70410066216957]],"pts":[{"x":121.48571681415929,"y":25.058614130434783,"z":95},{"x":121.48500884955752,"y":25.062228260869563,"z":146},{"x":121.47,"y":25.08125,"z":22},{"x":121.53598230088495,"y":25.02475543478261,"z":223},{"x":121.54582300884955,"y":25.074592391304346,"z":39},{"x":121.52118584070796,"y":25.057472826086958,"z":191},{"x":121.48083185840707,"y":25.02,"z":22},{"x":121.50653097345132,"y":25.089429347826087,"z":249},{"x":121.55,"y":25.09,"z":86},{"x":121.52614159292035,"y":25.067364130434783,"z":125}]}
*
*/
async function interp2Raster(ops, opt = {}) {
let dlng = 0.000979087
let dlat = 0.000906247
//check
if (!isearr(ops)) {
throw new Error(`ops`)
}
//keyX
let keyX = get(opt, 'keyX')
if (!isestr(keyX)) {
keyX = 'x'
}
//keyY
let keyY = get(opt, 'keyY')
if (!isestr(keyY)) {
keyY = 'y'
}
//keyZ
let keyZ = get(opt, 'keyZ')
if (!isestr(keyZ)) {
keyZ = 'z'
}
//dx
let dx = get(opt, 'dx')
if (!isnum(dx)) {
dx = dlng
}
dx = cdbl(dx)
//dy
let dy = get(opt, 'dy')
if (!isnum(dy)) {
dy = dlat
}
dy = cdbl(dy)
//dxAgr
let dxAgr = get(opt, 'dxAgr')
if (!isnum(dxAgr)) {
dxAgr = dy * 2
}
dxAgr = cdbl(dxAgr)
//dyAgr
let dyAgr = get(opt, 'dyAgr')
if (!isnum(dyAgr)) {
dyAgr = dy * 2
}
dyAgr = cdbl(dyAgr)
//modePick
let modePick = get(opt, 'modePick')
if (modePick !== 'min' && modePick !== 'max') {
modePick = 'min'
}
//funValid
let funValid = get(opt, 'funValid')
if (!isfun(funValid)) {
funValid = () => {
return true
}
}
//funKriging
let funKriging = get(opt, 'funKriging')
if (!isfun(funKriging)) {
funKriging = _kriging
}
//funAdjust
let funAdjust = get(opt, 'funAdjust')
if (!isfun(funAdjust)) {
funAdjust = (x, y, z) => {
return z
}
}
//xmin, xmax, ymin, ymax
let xmin = 1e20
let xmax = -1e20
let ymin = 1e20
let ymax = -1e20
each(ops, (m) => {
let v
v = get(m, keyX)
if (isNumber(v)) {
xmin = Math.min(v, xmin)
xmax = Math.max(v, xmax)
}
v = get(m, keyY)
if (isNumber(v)) {
ymin = Math.min(v, ymin)
ymax = Math.max(v, ymax)
}
})
// console.log('xmin', xmin, 'xmax', xmax)
// console.log('ymin', ymin, 'ymax', ymax)
//pts, 用dxAgr與dyAgr聚合提取點數據
let pts = aggregatePoints(ops, xmin, dxAgr, ymin, dyAgr, {
keyX,
keyY,
keyZ,
modePick,
})
// console.log('pts', take(pts, 5), size(pts))
//rg
let rg = await interp2Grid(pts, xmin, xmax, dx, ymin, ymax, dy, {
...opt, //須提供自定義參數給funKriging使用
funKriging,
funValid,
funAdjust,
returnGrid: true,
inverseKeyY: true,
})
// console.log('rg', rg)
// console.log('keys(rg)', keys(rg))
// keys(rg) [
// 'xnum', 'xmin', 'xmax',
// 'dx', 'ynum', 'ymin',
// 'ymax', 'dy', 'zmin',
// 'zmax', 'grds'
// ]
// console.log('rg.grds', rg.grds[0])
// rg.grds [
// 60.606247155177684, 62.92540214009434, 65.7312154911937, 68.9914947349822,
// 72.65919624422395, 76.67923534212574, 80.99475324860614, 85.55171796733839,
// ...
// ]
//save
rg.pts = pts
return rg
}
export default interp2Raster