From 97cb8f5735937ea9c98f03d87fea1d50213a82b8 Mon Sep 17 00:00:00 2001 From: thkim Date: Tue, 13 Jan 2026 09:31:54 +0900 Subject: [PATCH] =?UTF-8?q?fixed:=20=EB=88=84=EC=A0=81=EC=B9=A8=ED=95=98?= =?UTF-8?q?=EB=9F=89=EC=9D=B4=20=EB=B6=84=EB=AA=A8=EB=A1=9C=20=EA=B0=96?= =?UTF-8?q?=EB=8A=94=20=EC=97=B0=EC=82=B0=EC=97=90=EC=84=9C=20=20=EB=88=84?= =?UTF-8?q?=EC=A0=81=EC=B9=A8=ED=95=98=EB=9F=89=EC=9D=B4=200=EC=9D=BC=20?= =?UTF-8?q?=EB=95=8C=20=EB=B0=9C=EC=83=9D=ED=95=98=EB=8A=94=200=EC=9C=BC?= =?UTF-8?q?=EB=A1=9C=20=EB=82=98=EB=88=84=EA=B8=B0=20=EC=97=90=EB=9F=AC=20?= =?UTF-8?q?=EC=88=98=EC=A0=95?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- 20240119_backup/controller.py | 215 ----- .../settle_prediction_steps_main.py | 777 ------------------ ...ttle_prediction_steps_main.cpython-310.pyc | Bin 5055 -> 5969 bytes controller.py | 4 + settle_prediction_steps_main.py | 57 +- 5 files changed, 53 insertions(+), 1000 deletions(-) delete mode 100644 20240119_backup/controller.py delete mode 100644 20240119_backup/settle_prediction_steps_main.py diff --git a/20240119_backup/controller.py b/20240119_backup/controller.py deleted file mode 100644 index 0dc383b..0000000 --- a/20240119_backup/controller.py +++ /dev/null @@ -1,215 +0,0 @@ -""" -Title: Controller -Developer: -Sang Inn Woo, Ph.D. @ Incheon National University -Starting Date: 2022-11-10 -""" -import psycopg2 as pg2 -import sys -import numpy as np -import settle_prediction_steps_main -import matplotlib.pyplot as plt -import pdb; -''' -apptb_surset01 -cons_code: names of monitoring points - -apptb_surset02 -cons_code: names of monitoring points -amount_cum_sub: accumulated settlement -fill_height: height of surcharge fill -nod: number of date -''' - -def settlement_prediction(business_code, cons_code): - - # connect the database - #connection = pg2.connect("host=localhost dbname=sgis user=postgres password=postgres port=5434") # local - connection = pg2.connect("host=192.168.10.172 dbname=sgis_new user=sgis password=sgis port=5432") # ICTWay internal - - # set cursor - cursor = connection.cursor() - - # select monitoring data for the monitoring point - postgres_select_query = """SELECT (amount_cum_sub * -1), fill_height, nod FROM apptb_surset02 WHERE business_code='""" + business_code \ - + """' and cons_code='""" + cons_code + """' ORDER BY nod ASC""" - cursor.execute(postgres_select_query) - monitoring_record = cursor.fetchall() - # initialize time, surcharge, and settlement lists - time = [] - surcharge = [] - settlement = [] - - # fill lists - - for row in monitoring_record: - settlement.append(float(row[0])) - surcharge.append(float(row[1])) - time.append(float(row[2])) - - # convert lists to np arrays - settlement = np.array(settlement) - surcharge = np.array(surcharge) - time = np.array(time) - - # run the settlement prediction and get results - results = settle_prediction_steps_main.run_settle_prediction(point_name=cons_code, np_time=time, - np_surcharge=surcharge, np_settlement=settlement, - final_step_predict_percent=90, - additional_predict_percent=300, plot_show=False, - print_values=False, run_original_hyperbolic=True, - run_nonlinear_hyperbolic=True, - run_weighted_nonlinear_hyperbolic=True, - run_asaoka=True, run_step_prediction=True, - asaoka_interval=3) - - - # prediction method code - # 1: original hyperbolic method (쌍곡선법) - # 2: nonlinear hyperbolic method (비선형 쌍곡선법) - # 3: weighted nonlinear hyperbolic method (가중 비선형 쌍곡선법) - # 4: Asaoka method (아사오카법) - # 5: Step loading (단계성토 고려법) - - - ''' - time_hyper, sp_hyper_original, - time_hyper, sp_hyper_nonlinear, - time_hyper, sp_hyper_weight_nonlinear, - time_asaoka, sp_asaoka, - time[step_start_index[0]:], -sp_step[step_start_index[0]:], - ''' - - for i in range(5): - - # if there are prediction data for the given data point, delete it first - postgres_delete_query = """DELETE FROM apptb_pred02_no""" + str(i + 1) \ - + """ WHERE business_code='""" + business_code \ - + """' and cons_code='""" + cons_code + """'""" - cursor.execute(postgres_delete_query) - connection.commit() - - # get time and settlement arrays - time = results[2 * i] - predicted_settlement = results[2 * i + 1] - - # for each prediction time - for j in range(len(time)): - - # construct insert query - postgres_insert_query \ - = """INSERT INTO apptb_pred02_no""" + str(i + 1) + """ """ \ - + """(business_code, cons_code, prediction_progress_days, predicted_settlement, prediction_method) """ \ - + """VALUES (%s, %s, %s, %s, %s)""" - - # set data to insert - record_to_insert = (business_code, cons_code, time[j], predicted_settlement[j], i + 1) - - # execute the insert query - cursor.execute(postgres_insert_query, record_to_insert) - - # commit changes - connection.commit() - - -def read_database_and_plot(business_code, cons_code): - - # connect the database - # connection = pg2.connect("host=localhost dbname=postgres user=postgres password=lab36981 port=5432") - connection = pg2.connect("host=192.168.0.72 dbname=sgis user=sgis password=sgis port=5432") # ICTWay internal - - # set cursor - cursor = connection.cursor() - - # select monitoring data for the monitoring point - postgres_select_query = """SELECT * FROM apptb_surset02 WHERE business_code='""" + business_code \ - + """' and cons_code='""" + cons_code + """' ORDER BY nod ASC""" - cursor.execute(postgres_select_query) - monitoring_record = cursor.fetchall() - - # initialize time, surcharge, and settlement lists - time_monitored = [] - surcharge_monitored = [] - settlement_monitored = [] - - # fill lists - for row in monitoring_record: - time_monitored.append(float(row[2])) - settlement_monitored.append(float(row[6])) - surcharge_monitored.append(float(row[8])) - - # convert lists to np arrays - settlement_monitored = np.array(settlement_monitored) - surcharge_monitored = np.array(surcharge_monitored) - time_monitored = np.array(time_monitored) - - # prediction method code - # 0: original hyperbolic method - # 1: nonlinear hyperbolic method - # 2: weighted nonlinear hyperbolic method - # 3: Asaoka method - # 4: Step loading - # 5: temp - - # temporarily set the prediction method as 0 - postgres_select_query = """SELECT prediction_progress_days, predicted_settlement """ \ - + """FROM apptb_pred02_no""" + str(5) \ - + """ WHERE business_code='""" + business_code \ - + """' and cons_code='""" + cons_code \ - + """' ORDER BY prediction_progress_days ASC""" - - # select predicted data for the monitoring point - cursor.execute(postgres_select_query) - prediction_record = cursor.fetchall() - - # initialize time, surcharge, and settlement lists - time_predicted = [] - settlement_predicted = [] - - # fill lists - for row in prediction_record: - time_predicted.append(float(row[0])) - settlement_predicted.append(float(row[1])) - - # convert lists to np arrays - settlement_predicted = np.array(settlement_predicted) - time_predicted = np.array(time_predicted) - - # 그래프 크기, 서브 그래프 개수 및 비율 설정 - fig, axes = plt.subplots(2, 1, figsize=(8, 6), gridspec_kw={'height_ratios': [1, 3]}) - - # 성토고 그래프 표시 - axes[0].plot(time_monitored, surcharge_monitored, color='black', label='surcharge height') - - # 성토고 그래프 설정 - axes[0].set_ylabel("Surcharge height (m)", fontsize=10) - axes[0].set_xlim(left=0) - axes[0].set_xlim(right=np.max(time_predicted)) - axes[0].grid(color="gray", alpha=.5, linestyle='--') - axes[0].tick_params(direction='in') - - # 계측 및 예측 침하량 표시 - axes[1].scatter(time_monitored, -settlement_monitored, s=30, - facecolors='white', edgecolors='black', label='measured data') - - axes[1].plot(time_predicted, -settlement_predicted, - linestyle='--', color='red', label='Original Hyperbolic') - - axes[0].set_ylabel("Settlement (cm)", fontsize=10) - axes[1].set_xlim(left=0) - axes[1].set_xlim(right=np.max(time_predicted)) - - -# script to call: python3 controller.py [business_code] [cons_code] -# for example: python3 controller.py 221222SA0003 CONS001 -if __name__ == '__main__': - - args = sys.argv[1:] - business_code = args[0] - cons_code = args[1] - - settlement_prediction(business_code=business_code, cons_code=cons_code) - print("The settlement prediction is over.") - - #read_database_and_plot(business_code=business_code, cons_code=cons_code) - #print("Visualization is over.") \ No newline at end of file diff --git a/20240119_backup/settle_prediction_steps_main.py b/20240119_backup/settle_prediction_steps_main.py deleted file mode 100644 index 1435110..0000000 --- a/20240119_backup/settle_prediction_steps_main.py +++ /dev/null @@ -1,777 +0,0 @@ -""" -Title: Soft ground settlement prediction -Developer: -Sang Inn Woo, Ph.D. @ Incheon National University -Kwak Taeyoung, Ph.D. @ KICT - -Starting Date: 2022-08-11 -Abstract: -This main objective of this code is to predict -time vs. (consolidation) settlement curves of soft clay ground. -""" - -# ================= -# Import 섹션 -# ================= -import os.path -import numpy as np -import pandas as pd -import matplotlib.pyplot as plt -from scipy.optimize import least_squares -from scipy.interpolate import interp1d - - -# ================= -# Function 섹션 -# ================= - -# 주어진 계수를 이용하여 쌍곡선 시간-침하 곡선 반환 -def generate_data_hyper(px, pt): - return pt / (px[0] * pt + px[1]) - -# 주어진 계수를 이용하여 아사오카 시간-침하 곡선 반환 -def generate_data_asaoka(px, pt, dt): - return (px[1] / (1 - px[0])) * (1 - (px[0] ** (pt / dt))) - -# 회귀식과 측정치와의 잔차 반환 (비선형 쌍곡선) -def fun_hyper_nonlinear(px, pt, py): - return pt / (px[0] * pt + px[1]) - py - -# 회귀식과 측정치와의 잔차 반환 (가중 비선형 쌍곡선) -def fun_hyper_weight_nonlinear(px, pt, py, pw): - return (pt / (px[0] * pt + px[1]) - py) * pw - -# 회귀식과 측정치와의 잔차 반환 (기존 쌍곡선) -def fun_hyper_original(px, pt, py): - return px[0] * pt + px[1] - pt / py - -# 회귀식과 측정치와의 잔차 반환 (아사오카) -def fun_asaoka(px, ps_b, ps_a): - return px[0] * ps_b + px[1] - ps_a - - -# RMSE 산정 -def fun_rmse(py1, py2): - mse = np.square(np.subtract(py1, py2)).mean() - return np.sqrt(mse) - - -def run_settle_prediction_from_file(input_file, output_dir, - final_step_predict_percent, additional_predict_percent, - plot_show, - print_values, - run_original_hyperbolic='True', - run_nonlinear_hyperbolic='True', - run_weighted_nonlinear_hyperbolic='True', - run_asaoka='True', - run_step_prediction='True', - asaoka_interval=3): - - # 현재 파일 이름 출력 - print("Working on " + input_file) - - # CSV 파일 읽기 - data = pd.read_csv(input_file, encoding='euc-kr') - - # 시간, 침하량, 성토고 배열 생성 - time = data['Time'].to_numpy() - settle = data['Settlement'].to_numpy() - surcharge = data['Surcharge'].to_numpy() - - run_settle_prediction(point_name=input_file, np_time=time, np_surcharge=surcharge, np_settlement=settle, - final_step_predict_percent=final_step_predict_percent, - additional_predict_percent=additional_predict_percent, plot_show=plot_show, - print_values=print_values, - run_original_hyperbolic=run_original_hyperbolic, - run_nonlinear_hyperbolic=run_nonlinear_hyperbolic, - run_weighted_nonlinear_hyperbolic=run_weighted_nonlinear_hyperbolic, - run_asaoka=run_asaoka, - run_step_prediction=run_step_prediction, - asaoka_interval=asaoka_interval) - -def run_settle_prediction(point_name, - np_time, np_surcharge, np_settlement, - final_step_predict_percent, additional_predict_percent, - plot_show, - print_values, - run_original_hyperbolic='True', - run_nonlinear_hyperbolic='True', - run_weighted_nonlinear_hyperbolic='True', - run_asaoka = 'True', - run_step_prediction='True', - asaoka_interval = 5): - - # ==================== - # 파일 읽기, 데이터 설정 - # ==================== - - # 시간, 침하량, 성토고 배열 생성 - time = np_time - settle = np_settlement - surcharge = np_surcharge - - # 마지막 계측 데이터 index + 1 파악 - final_index = time.size - - # ================= - # 성토 단계 구분 - # ================= - - # 성토 단계 시작 index 리스트 초기화 - step_start_index = [0] - - # 성토 단계 끝 index 리스트 초기화 - step_end_index = [] - - # 현재 성토고 설정 - current_surcharge = surcharge[0] - - # 단계 시작 시점 초기화 - step_start_date = 0 - # 모든 시간-성토고 데이터에서 순차적으로 확인 - for index in range(len(surcharge)): - - # 만일 성토고의 변화가 있을 경우, - if surcharge[index] > current_surcharge*1.05 or surcharge[index] < current_surcharge*0.95: - step_end_index.append(index) - step_start_index.append(index) - current_surcharge = surcharge[index] - - # 마지막 성토 단계 끝 index 추가 - step_end_index.append(len(surcharge) - 1) - - # ================= - # 성토 단계 조정 - # ================= - # 성토고 유지 기간이 매우 짧을 경우, 해석 단계에서 제외 - - # 조정 성토 시작 및 끝 인덱스 리스트 초기화 - step_start_index_adjust = [] - step_end_index_adjust = [] - - # 각 성토 단계 별로 분석 - for i in range(0, len(step_start_index)): - - # 현 단계 성토 시작일 / 끝일 파악 - step_start_date = time[step_start_index[i]] - step_end_date = time[step_end_index[i]] - - # 현 성토고 유지 일수 및 데이터 개수 파악 - step_span = step_end_date - step_start_date - step_data_num = step_end_index[i] - step_start_index[i] + 1 - - # 성토고 유지일 및 데이터 개수 기준 적용 - if step_span > 30 and step_data_num > 5: - step_start_index_adjust.append((step_start_index[i])) - step_end_index_adjust.append((step_end_index[i])) - - # 성토 시작 및 끝 인덱스 리스트 업데이트 - step_start_index = step_start_index_adjust - step_end_index = step_end_index_adjust - - # 침하 예측을 수행할 단계 설정 (현재 끝에서 2단계 이용) - step_start_index = step_start_index[-2:] - step_end_index = step_end_index[-2:] - - # 성토 단계 횟수 파악 및 저장 - num_steps = len(step_start_index) - - # =========================== - # 최종 단계 데이터 사용 범위 조정 - # =========================== - - # 데이터 사용 퍼센트에 해당하는 기간 계산 - final_step_end_date = time[-1] - final_step_start_date = time[step_start_index[num_steps - 1]] - final_step_period = final_step_end_date - final_step_start_date - final_step_predict_end_date = final_step_start_date + final_step_period * final_step_predict_percent / 100 - - # 데이터 사용 끝 시점 인덱스 초기화 - final_step_predict_end_index = -1 - - # 데이터 사용 끝 시점 인덱스 검색 - count = 0 - for day in time: - count = count + 1 - if day > final_step_predict_end_date: - final_step_predict_end_index = count - 1 - break - - # 마지막 성토 단계, 마지막 계측 시점 인덱스 업데이트 - final_step_monitor_end_index = step_end_index[num_steps - 1] - step_end_index[num_steps - 1] = final_step_predict_end_index - - # ================= - # 추가 예측 구간 반영 - # ================= - - # 추가 예측 일 입력 (현재 전체 계측일 * 계수) - add_days = (additional_predict_percent / 100) * time[-1] - - # 마지막 성토고 및 마지막 계측일 저장 - final_surcharge = surcharge[final_index - 1] - final_time = time[final_index - 1] - - # 추가 시간 및 성토고 배열 설정 (100개의 시점 설정) - time_add = np.linspace(final_time + 1, final_time + add_days, 100) - surcharge_add = np.ones(100) * final_surcharge - - # 기존 시간 및 성토고 배열에 붙이기 - time = np.append(time, time_add) - surcharge = np.append(surcharge, surcharge_add) - - # 마지막 인덱스값 재조정 - final_index = time.size - - - - - - - - - - # ========================================== - # Settlement Prediction (Step + Hyperbolic) - # ========================================== - - # 예측 침하량 초기화 - sp_step = np.zeros(time.size) - - # 각 단계별로 진행 - for i in range(0, num_steps): - - # 각 단계별 계측 시점과 계측 침하량 배열 생성 - tm_this_step = time[step_start_index[i]:step_end_index[i]] - sm_this_step = settle[step_start_index[i]:step_end_index[i]] - - # 이전 단계까지 예측 침하량 중 현재 단계에 해당하는 부분 추출 - sp_this_step = sp_step[step_start_index[i]:step_end_index[i]] - - # 현재 단계 시작 부터 끝까지 시간 데이터 추출 - tm_to_end = time[step_start_index[i]:final_index] - - # 기존 예측 침하량에 대한 보정 - sm_this_step = sm_this_step - sp_this_step - - # 초기 시점 및 침하량 산정 - t0_this_step = tm_this_step[0] - s0_this_step = sm_this_step[0] - - # 초기 시점에 대한 시간 조정 - tm_this_step = tm_this_step - t0_this_step - tm_to_end = tm_to_end - t0_this_step - - # 초기 침하량에 대한 침하량 조정 - sm_this_step = sm_this_step - s0_this_step - - # 침하 곡선 계수 초기화 - x0 = np.ones(2) - - # 회귀분석 시행 - res_lsq_hyper_nonlinear \ - = least_squares(fun_hyper_nonlinear, x0, - args=(tm_this_step, sm_this_step)) - - # 쌍곡선 계수 저장 및 출력 - x_step = res_lsq_hyper_nonlinear.x - if print_values: - print(x_step) - - # 현재 단계 예측 침하량 산정 (침하 예측 끝까지) - sp_to_end_update = generate_data_hyper(x_step, tm_to_end) - - # 예측 침하량 업데이트 - sp_step[step_start_index[i]:final_index] = \ - sp_step[step_start_index[i]:final_index] + sp_to_end_update + s0_this_step - - - - - ''' - # ====================================== - # Settlement Prediction (Step + Asaoka) - # ====================================== - - # TODO: Modify this - - # 예측 침하량 초기화 - sp_step_asaoka = np.zeros(time.size) - - # 각 단계별로 진행 - for i in range(0, num_steps): - - # 각 단계별 계측 시점과 계측 침하량 배열 생성 - tm_this_step = time[step_start_index[i]:step_end_index[i]] - sm_this_step = settle[step_start_index[i]:step_end_index[i]] - - # 이전 단계 까지 예측 침하량 중 현재 단계에 해당 하는 부분 추출 - sp_this_step = sp_step[step_start_index[i]:step_end_index[i]] - - # 현재 단계 시작 부터 끝까지 시간 데이터 추출 - tm_to_end = time[step_start_index[i]:final_index] - - # 기존 예측 침하량에 대한 보정 - sm_this_step = sm_this_step - sp_this_step - - # 초기 시점 및 침하량 산정 - t0_this_step = tm_this_step[0] - s0_this_step = sm_this_step[0] - - # 초기 시점에 대한 시간 조정 - tm_this_step = tm_this_step - t0_this_step - tm_to_end = tm_to_end - t0_this_step - - # 초기 침하량에 대한 침하량 조정 - sm_this_step = sm_this_step - s0_this_step - - - - - - - - # 등간격 데이터 생성을 위한 Interpolation 함수 설정 - inter_fn = interp1d(tm_this_step, sm_this_step, kind='cubic') - - # 데이터 구축 간격 및 그에 해당하는 데이터 포인트 개수 설정 - num_data = int(tm_this_step[-1] / asaoka_interval) - - # 등간격 시간 및 침하량 데이터 설정 - tm_this_step_inter = np.linspace(0, tm_this_step[-1], num=num_data, endpoint=True) - sm_this_step_inter = inter_fn(tm_this_step_inter) - - # 이전 이후 등간격 침하량 배열 구축 - sm_this_step_before = sm_this_step_inter[0:-2] - sm_this_step_after = sm_this_step_inter[1:-1] - - - - - # Least square 변수 초기화 - x0 = np.ones(2) - - # Least square 분석을 통한 침하 곡선 계수 결정 - res_lsq_asaoka = least_squares(fun_asaoka, x0, args=(sm_this_step_before, sm_this_step_after)) - - # 기존 쌍곡선 법 계수 저장 및 출력 - x_step_asaoka = res_lsq_asaoka.x - if print_values: - print(x_step_asaoka) - - # 현재 단계 예측 침하량 산정 (침하 예측 끝까지) - sp_to_end_update = generate_data_asaoka(x_step_asaoka, tm_to_end, asaoka_interval) - - - # 예측 침하량 업데이트 - sp_step_asaoka[step_start_index[i]:final_index] = \ - sp_step_asaoka[step_start_index[i]:final_index] + sp_to_end_update + s0_this_step - - ''' - - # ========================================================= - # Settlement prediction (nonliner, weighted nonlinear and original hyperbolic) - # ========================================================= - - # 성토 마지막 데이터 추출 - tm_hyper = time[step_start_index[num_steps - 1]:step_end_index[num_steps - 1]] - sm_hyper = settle[step_start_index[num_steps - 1]:step_end_index[num_steps - 1]] - - # 현재 단계 시작 부터 끝까지 시간 데이터 추출 - time_hyper = time[step_start_index[num_steps - 1]:final_index] - - # 초기 시점 및 침하량 산정 - t0_hyper = tm_hyper[0] - s0_hyper = sm_hyper[0] - - # 초기 시점에 대한 시간 조정 - tm_hyper = tm_hyper - t0_hyper - time_hyper = time_hyper - t0_hyper - - # 초기 침하량에 대한 침하량 조정 - sm_hyper = sm_hyper - s0_hyper - - # 회귀분석 시행 (비선형 쌍곡선) - x0 = np.ones(2) - res_lsq_hyper_nonlinear = least_squares(fun_hyper_nonlinear, x0, - args=(tm_hyper, sm_hyper)) - # 비선형 쌍곡선 법 계수 저장 및 출력 - x_hyper_nonlinear = res_lsq_hyper_nonlinear.x - if print_values: - print(x_hyper_nonlinear) - - # 가중 비선형 쌍곡선 가중치 산정 - weight = tm_hyper / np.sum(tm_hyper) - - # 회귀분석 시행 (가중 비선형 쌍곡선) - x0 = np.ones(2) - res_lsq_hyper_weight_nonlinear = least_squares(fun_hyper_weight_nonlinear, x0, - args=(tm_hyper, sm_hyper, weight)) - # 비선형 쌍곡선 법 계수 저장 및 출력 - x_hyper_weight_nonlinear = res_lsq_hyper_weight_nonlinear.x - if print_values: - print(x_hyper_weight_nonlinear) - - # 회귀분석 시행 (기존 쌍곡선법) - (0, 0)에 해당하는 초기 데이터를 제외하고 회귀분석 실시 - x0 = np.ones(2) - res_lsq_hyper_original = least_squares(fun_hyper_original, x0, - args=(tm_hyper[1:], sm_hyper[1:])) - # 기존 쌍곡선 법 계수 저장 및 출력 - x_hyper_original = res_lsq_hyper_original.x - if print_values: - print(x_hyper_original) - - # 현재 단계 예측 침하량 산정 (침하 예측 끝까지) - sp_hyper_nonlinear = generate_data_hyper(x_hyper_nonlinear, time_hyper) - sp_hyper_weight_nonlinear = generate_data_hyper(x_hyper_weight_nonlinear, time_hyper) - sp_hyper_original = generate_data_hyper(x_hyper_original, time_hyper) - - # 예측 침하량 산정 - sp_hyper_nonlinear = sp_hyper_nonlinear + s0_hyper - sp_hyper_weight_nonlinear = sp_hyper_weight_nonlinear + s0_hyper - sp_hyper_original = sp_hyper_original + s0_hyper - time_hyper = time_hyper + t0_hyper - - - - - - - - # =============================== - # Settlement prediction (Asaoka) - # =============================== - - # 성토 마지막 데이터 추출 - tm_asaoka = time[step_start_index[num_steps - 1]:step_end_index[num_steps - 1]] - sm_asaoka = settle[step_start_index[num_steps - 1]:step_end_index[num_steps - 1]] - - # 현재 단계 시작 부터 끝까지 시간 데이터 추출 - time_asaoka = time[step_start_index[num_steps - 1]:final_index] - - # 초기 시점 및 침하량 산정 - t0_asaoka = tm_asaoka[0] - s0_asaoka = sm_asaoka[0] - - # 초기 시점에 대한 시간 조정 - tm_asaoka = tm_asaoka - t0_asaoka - time_asaoka = time_asaoka - t0_asaoka - - # 초기 침하량에 대한 침하량 조정 - sm_asaoka = sm_asaoka - s0_asaoka - - # 등간격 데이터 생성을 위한 Interpolation 함수 설정 - inter_fn = interp1d(tm_asaoka, sm_asaoka, kind='cubic') - - # 데이터 구축 간격 및 그에 해당하는 데이터 포인트 개수 설정 - num_data = int(tm_asaoka[-1] / asaoka_interval) - - # 등간격 시간 및 침하량 데이터 설정 - tm_asaoka_inter = np.linspace(0, tm_asaoka[-1], num=num_data, endpoint=True) - sm_asaoka_inter = inter_fn(tm_asaoka_inter) - - # 이전 이후 등간격 침하량 배열 구축 - sm_asaoka_before = sm_asaoka_inter[0:-2] - sm_asaoka_after = sm_asaoka_inter[1:-1] - - # Least square 변수 초기화 - x0 = np.ones(2) - - # Least square 분석을 통한 침하 곡선 계수 결정 - res_lsq_asaoka = least_squares(fun_asaoka, x0, args=(sm_asaoka_before, sm_asaoka_after)) - - # 기존 쌍곡선 법 계수 저장 및 출력 - x_asaoka = res_lsq_asaoka.x - if print_values: - print(x_asaoka) - - # 현재 단계 예측 침하량 산정 (침하 예측 끝까지) - sp_asaoka = generate_data_asaoka(x_asaoka, time_asaoka, asaoka_interval) - - # 예측 침하량 산정 - sp_asaoka = sp_asaoka + s0_asaoka - time_asaoka = time_asaoka + t0_asaoka - - - - - - - - - - - - - - - - # ============================== - # Post-Processing #1 : 에러 산정 - # ============================== - - # RMSE 계산 데이터 구간 설정 (계측) - sm_rmse = settle[final_step_predict_end_index:final_step_monitor_end_index] - - # RMSE 계산 데이터 구간 설정 (단계) - sp_step_rmse = sp_step[final_step_predict_end_index:final_step_monitor_end_index] - - # RMSE 계산 데이터 구간 설정 (쌍곡선) - sp_hyper_nonlinear_rmse = sp_hyper_nonlinear[final_step_predict_end_index - step_start_index[num_steps - 1]: - final_step_predict_end_index - step_start_index[num_steps - 1] + - final_step_monitor_end_index - final_step_predict_end_index] - sp_hyper_weight_nonlinear_rmse \ - = sp_hyper_weight_nonlinear[final_step_predict_end_index - step_start_index[num_steps - 1]: - final_step_predict_end_index - step_start_index[num_steps - 1] + - final_step_monitor_end_index - final_step_predict_end_index] - sp_hyper_original_rmse = sp_hyper_original[final_step_predict_end_index - step_start_index[num_steps - 1]: - final_step_predict_end_index - step_start_index[num_steps - 1] + - final_step_monitor_end_index - final_step_predict_end_index] - - # RMSE 계산 데이터 구간 설정 (아사오카) - sp_asaoka_rmse = sp_asaoka[final_step_predict_end_index - step_start_index[num_steps - 1]: - final_step_predict_end_index - step_start_index[num_steps - 1] + - final_step_monitor_end_index - final_step_predict_end_index] - - # RMSE 산정 (단계, 비선형 쌍곡선, 기존 쌍곡선) - rmse_step = fun_rmse(sm_rmse, sp_step_rmse) - rmse_hyper_nonlinear = fun_rmse(sm_rmse, sp_hyper_nonlinear_rmse) - rmse_hyper_weight_nonlinear = fun_rmse(sm_rmse, sp_hyper_weight_nonlinear_rmse) - rmse_hyper_original = fun_rmse(sm_rmse, sp_hyper_original_rmse) - rmse_asaoka = fun_rmse(sm_rmse, sp_asaoka_rmse) - - # RMSE 출력 (단계, 비선형 쌍곡선, 기존 쌍곡선) - if print_values: - print("RMSE (Nonlinear Hyper + Step): %0.3f" % rmse_step) - print("RMSE (Nonlinear Hyperbolic): %0.3f" % rmse_hyper_nonlinear) - print("RMSE (Weighted Nonlinear Hyperbolic): %0.3f" % rmse_hyper_weight_nonlinear) - print("RMSE (Original Hyperbolic): %0.3f" % rmse_hyper_original) - print("RMSE (Asaoka): %0.3f" % rmse_asaoka) - - # (최종 계측 침하량 - 예측 침하량) 계산 - final_error_step = np.abs(settle[-1] - sp_step_rmse[-1]) - final_error_hyper_nonlinear = np.abs(settle[-1] - sp_hyper_nonlinear_rmse[-1]) - final_error_hyper_weight_nonlinear = np.abs(settle[-1] - sp_hyper_weight_nonlinear_rmse[-1]) - final_error_hyper_original = np.abs(settle[-1] - sp_hyper_original_rmse[-1]) - final_error_asaoka = np.abs(settle[-1] - sp_asaoka_rmse[-1]) - - # (최종 계측 침하량 - 예측 침하량) 출력 (단계, 비선형 쌍곡선, 기존 쌍곡선) - if print_values: - print("Error in Final Settlement (Nonlinear Hyper + Step): %0.3f" % final_error_step) - print("Error in Final Settlement (Nonlinear Hyperbolic): %0.3f" % final_error_hyper_nonlinear) - print("Error in Final Settlement (Weighted Nonlinear Hyperbolic): %0.3f" % final_error_hyper_weight_nonlinear) - print("Error in Final Settlement (Original Hyperbolic): %0.3f" % final_error_hyper_original) - print("Error in Final Settlement (Asaoka): %0.3f" % final_error_asaoka) - - - - - # ========================================== - # Post-Processing #2 : 그래프 작성 - # ========================================== - - # 만약 그래프 도시가 필요할 경우, - if plot_show: - - # 그래프 크기, 서브 그래프 개수 및 비율 설정 - fig, axes = plt.subplots(2, 1, figsize=(12, 9), gridspec_kw={'height_ratios': [1, 3]}) - - # 성토고 그래프 표시 - axes[0].plot(time, surcharge, color='black', label='surcharge height') - - # 성토고 그래프 설정 - axes[0].set_ylabel("Surcharge height (m)", fontsize=15) - axes[0].set_xlim(left=0) - axes[0].grid(color="gray", alpha=.5, linestyle='--') - axes[0].tick_params(direction='in') - - # 계측 및 예측 침하량 표시 - axes[1].scatter(time[0:settle.size], -settle, s=50, - facecolors='white', edgecolors='black', label='measured data') - axes[1].plot(time_hyper, -sp_hyper_original, - linestyle='--', color='red', label='Original Hyperbolic') - axes[1].plot(time_hyper, -sp_hyper_nonlinear, - linestyle='--', color='green', label='Nonlinear Hyperbolic') - axes[1].plot(time_hyper, -sp_hyper_weight_nonlinear, - linestyle='--', color='blue', label='Nonlinear Hyperbolic (Weighted)') - axes[1].plot(time_asaoka, -sp_asaoka, - linestyle='--', color='orange', label='Asaoka') - axes[1].plot(time[step_start_index[0]:], -sp_step[step_start_index[0]:], - linestyle='--', color='navy', label='Nonlinear + Step Loading') - - # 침하량 그래프 설정 - axes[1].set_xlabel("Time (day)", fontsize=15) - axes[1].set_ylabel("Settlement (cm)", fontsize=15) - axes[1].set_ylim(top=0) - axes[1].set_ylim(bottom=-1.5 * settle.max()) - axes[1].set_xlim(left=0) - axes[1].grid(color="gray", alpha=.5, linestyle='--') - axes[1].tick_params(direction='in') - - # 범례 표시 - axes[1].legend(loc=1, ncol=3, frameon=True, fontsize=10) - - # 예측 데이터 사용 범위 음영 처리 - 단계성토 - plt.axvspan(time[step_start_index[0]], final_step_predict_end_date, - alpha=0.1, color='grey', hatch='//') - - # 예측 데이터 사용 범위 음영 처리 - 기존 및 비선형 쌍곡선 - plt.axvspan(final_step_start_date, final_step_predict_end_date, - alpha=0.1, color='grey', hatch='\\') - - # 예측 데이터 사용 범위 표시 화살표 세로 위치 설정 - arrow1_y_loc = 1.3 * min(-settle) - arrow2_y_loc = 1.4 * min(-settle) - - # 화살표 크기 설정 - arrow_head_width = 0.03 * max(settle) - arrow_head_length = 0.01 * max(time) - - # 예측 데이터 사용 범위 화살표 처리 - 단계성토 - axes[1].arrow(time[step_start_index[0]], arrow1_y_loc, - final_step_predict_end_date - time[step_start_index[0]], 0, - head_width=arrow_head_width, head_length=arrow_head_length, - color='black', length_includes_head='True') - axes[1].arrow(final_step_predict_end_date, arrow1_y_loc, - time[step_start_index[0]] - final_step_predict_end_date, 0, - head_width=arrow_head_width, head_length=arrow_head_length, - color='black', length_includes_head='True') - - # 예측 데이터 사용 범위 화살표 처리 - 기존 및 비선형 쌍곡선 - axes[1].arrow(final_step_start_date, arrow2_y_loc, - final_step_predict_end_date - final_step_start_date, 0, - head_width=arrow_head_width, head_length=arrow_head_length, - color='black', length_includes_head='True') - axes[1].arrow(final_step_predict_end_date, arrow2_y_loc, - final_step_start_date - final_step_predict_end_date, 0, - head_width=arrow_head_width, head_length=arrow_head_length, - color='black', length_includes_head='True') - - # Annotation 표시용 공간 설정 - space = max(time) * 0.01 - - # 예측 데이터 사용 범위 범례 표시 - 단계성토 - plt.annotate('Data Range Used (Nonlinear + Step Loading)', xy=(final_step_predict_end_date, arrow1_y_loc), - xytext=(final_step_predict_end_date + space, arrow1_y_loc), - horizontalalignment='left', verticalalignment='center') - - # 예측 데이터 사용 범위 범례 표시 - 기존 및 비선형 쌍곡선 - plt.annotate('Data Range Used (Hyperbolic and Asaoka)', xy=(final_step_predict_end_date, arrow1_y_loc), - xytext=(final_step_predict_end_date + space, arrow2_y_loc), - horizontalalignment='left', verticalalignment='center') - - # RMSE 산정 범위 표시 화살표 세로 위치 설정 - arrow3_y_loc = 0.55 * min(-settle) - - # RMSE 산정 범위 화살표 표시 - axes[1].arrow(final_step_predict_end_date, arrow3_y_loc, - final_step_end_date - final_step_predict_end_date, 0, - head_width=arrow_head_width, head_length=arrow_head_length, - color='black', length_includes_head='True') - axes[1].arrow(final_step_end_date, arrow3_y_loc, - final_step_predict_end_date - final_step_end_date, 0, - head_width=arrow_head_width, head_length=arrow_head_length, - color='black', length_includes_head='True') - - # RMSE 산정 범위 세로선 설정 - axes[1].axvline(x=final_step_end_date, color='silver', linestyle=':') - - # RMSE 산정 범위 범례 표시 - plt.annotate('RMSE Estimation Section', xy=(final_step_end_date, arrow3_y_loc), - xytext=(final_step_end_date + space, arrow3_y_loc), - horizontalalignment='left', verticalalignment='center') - - # RMSE 출력 - mybox = {'facecolor': 'white', 'edgecolor': 'black', 'boxstyle': 'round', 'alpha': 0.2} - plt.text(max(time) * 1.04, 0.20 * min(-settle), - r"$\bf{Root\ Mean\ Squared\ Error}$" - + "\n" + "Original Hyperbolic: %0.3f" % rmse_hyper_original - + "\n" + "Nonlinear Hyperbolic: %0.3f" % rmse_hyper_nonlinear - + "\n" + "Nonlinear Hyperbolic (Weighted): %0.3f" % rmse_hyper_weight_nonlinear - + "\n" + "Asaoka: %0.3f" % rmse_asaoka - + "\n" + "Nonlinear + Step Loading: %0.3f" % rmse_step, - color='r', horizontalalignment='right', - verticalalignment='top', fontsize='10', bbox=mybox) - - # (최종 계측 침하량 - 예측값) 출력 - plt.text(max(time) * 1.04, 0.55 * min(-settle), - r"$\bf{Error\ in\ Final\ Settlement}$" - + "\n" + "Original Hyperbolic: %0.3f" % final_error_hyper_original - + "\n" + "Nonlinear Hyperbolic: %0.3f" % final_error_hyper_nonlinear - + "\n" + "Nonlinear Hyperbolic (Weighted): %0.3f" % final_error_hyper_weight_nonlinear - + "\n" + "Asaoka: %0.3f" % final_error_asaoka - + "\n" + "Nonlinear + Step Loading: %0.3f" % final_error_step, - color='r', horizontalalignment='right', - verticalalignment='top', fontsize='10', bbox=mybox) - - # 파일 이름만 추출 - filename = os.path.basename(point_name) - - # 그래프 제목 표시 - plt.title(filename + ": up to %i%% data used in the final step" % final_step_predict_percent) - - # 그래프 저장 (SVG 및 PNG) - # plt.savefig(output_dir + '/' + filename +' %i percent (SVG).svg' %final_step_predict_percent, bbox_inches='tight') - #plt.savefig(output_dir + '/' + filename + ' %i percent (PNG).png' % final_step_predict_percent, bbox_inches='tight') - - # 그래프 출력 - if plot_show: - plt.show() - - # 그래프 닫기 (메모리 소모 방지) - plt.close() - - # 예측 완료 표시 - print("Settlement prediction is done for " + filename + - " with " + str(final_step_predict_percent) + "% data usage") - - # 단계 성토 고려 여부 표시 - is_multi_step = True - if len(step_start_index) == 1: - is_multi_step = False - - # 반환 - - axes[1].plot(time_hyper, -sp_hyper_original, - linestyle='--', color='red', label='Original Hyperbolic') - axes[1].plot(time_hyper, -sp_hyper_nonlinear, - linestyle='--', color='green', label='Nonlinear Hyperbolic') - axes[1].plot(time_hyper, -sp_hyper_weight_nonlinear, - linestyle='--', color='blue', label='Nonlinear Hyperbolic (Weighted)') - axes[1].plot(time_asaoka, -sp_asaoka, - linestyle='--', color='orange', label='Asaoka') - axes[1].plot(time[step_start_index[0]:], -sp_step[step_start_index[0]:], - linestyle='--', color='navy', label='Nonlinear + Step Loading') - - return [time_hyper, sp_hyper_original, - time_hyper, sp_hyper_nonlinear, - time_hyper, sp_hyper_weight_nonlinear, - time_asaoka, sp_asaoka, - time[step_start_index[0]:], sp_step[step_start_index[0]:], - rmse_hyper_original, - rmse_hyper_nonlinear, - rmse_hyper_weight_nonlinear, - rmse_asaoka, - rmse_step, - final_error_hyper_original, - final_error_hyper_nonlinear, - final_error_hyper_weight_nonlinear, - final_error_asaoka, - final_error_step] - - -''' -run_settle_prediction(input_file='data/2-5_No.39.csv', - output_dir='output', - final_step_predict_percent=50, - additional_predict_percent=100, - plot_show=True, - print_values=True, - run_original_hyperbolic=True, - run_nonlinear_hyperbolic=True, - run_weighted_nonlinear_hyperbolic=True, - run_asaoka=True, - run_step_prediction=True, - asaoka_interval=3, - settle_unit='cm') -''' \ No newline at end of file diff --git a/__pycache__/settle_prediction_steps_main.cpython-310.pyc b/__pycache__/settle_prediction_steps_main.cpython-310.pyc index ce6fa54a79f8fa5febcf850c4e7e06818118e39d..0c2154ed452bbe941470e70e28c5191cfa395c3e 100644 GIT binary patch delta 2236 zcmZ`)O>7%Q6yBNr$;NgZ$97{UL1Ie|jVa`ZQZ-52{QtFSenbTosVcL{Y|?F6s3^IS3lb`v5aPmtOC<^th((Aa>Ioqs4i%!jnT_*<%37XzGxOft z?|tvh+pkVM?2H#18X^)rO<%uZe?RzX{4Z2hV5!1?6#i}}3``=ZMPew;CEOxqqzd+= zZNxA^tC+|P@5mLIAttv;m!(~0mLjlv#E0L$`0~e}eq}$O-^BrTN`bm@-J(`b32^Ba z{JlUWmQ!2{B4WzB=CY;iN;K}z;Kiz|frm6oA|25Lh)XgIwTC0zM>z!xft7Spm7*-D&7sQo9053jOPc~AZQ$oTfIV~|1*zqACx3mjWYN4+bPBd(nY0TI_8E^Mn)ZK@ zhfnRwVKzfb1P{sq4>Scv3%Im#_yzt2<=N_K?hV6WAleMmk+6?c*eh^)P)$3Cwt0{- zhKNMucc9unli28@l4xsjj~0H#!a)ZZXGvDp5dBnrGnLOJQ779($A%wLUot|*>pAo zwM1>lqxIqbep!Tc)Ri{#30}!ZNO_=SQhlKBaxZb*3}yWcqOJb#O8?j~*b>S)S#Z0c zl-&HH7|p{@KrM*y)S0eXTs4QeJI-+yNb^_K2~|M!KUce_|C4Z?ThTSwn&38^|7H|H z)Bu5IU|8Y>RD0$@;KeN(z`NT0ZmmwKDx`ryE=b(NtZNH_)aa?(3D(UsvIFaTn#fXZ z5AC44nC3i%fR50D9tJ`%X zx`Qew%W<#XOR|%!pPgOCBz8xxoMY#1L@Vd7NKUsRd8xa2a~(;Q0VI{v@cJe;aOkzd z3rVIWuK&73K<~`~CmqfkXhQ;vln69^Pu>fN+l->nD14O{V=uFwEnFD{saHrtAe|S| zFdJm2phpg}Au$iJVKEQ;??>9@5vKbeMOvz(td2!lGf}b<^wtQk^vW0;5p0aJQ87=j z5jGC}Y=n)8d5n#+30RMMX*d}r$4E1TKUpt%$H1~4{Ag{`7OG5n>J1ILTevqxG9-IP zsa#-_>;jx#WETlqB)WLzG%OTSFXT&*&^?xc<7VJ^nqAuCc!tdY$Fp9B%>o#0ORmhZ zIW_@d&T()skFg#IvbKK3W9`pn1BJ&IsMK^;tUxYUPN12klI0NE0#}Uc z{>L$+J`G}kZqbDQFqS^e`z6itfj1c1$`2#2@dhe2-flYIp)nqklVfqEXXhx~40-uB%v!KlvB0>JCHz delta 1267 zcmZWpOHUI~6ux&lopzv|7N(R!@fBjEX;d1GAT6|%mslh*?nn$NPQfa+&7CoVxf2#d zjWMq7#wEHjtc>voxNzf!L^m$xAFweoCMLQN&$&<@GR-~byN`3e^W8bquNOX^GwR`R z3xTI{BU61l@ZR`Erxlp0Fveim2{4>;V$m+}F#<^GBAp~jZjnvjETgc$jb@Xi6i;ca zz;Kub3rrKwl!djQm@YvXvfS~3X+DVkFeFVMfL#&`$$%mO6io~TC`n?%S`27YK#PHr zTT-QbEPx9^FrH#H`Q;g1Z1}p;X8_w+NgrilaK|dy7u+LWsVDkw5r!)=`{t4CjnzGM zxUGwuY9f3Vpp@A|SyVh$N7DzvDfOT<7J-9vO5@Ry4BypLYbn@}<(;Oi9q_U8@UrAZ z^#C~jHg+#40he>Y%nUesmL^(i15)Jb#ztlh0V1B zAHn8#M|ocQA|K(S3F40Mg7gKR=S7(3oi^CXTVcz98pd?uv_qs0^=V16PTg_Ge@}ys z-d!8FB3AUJ>Q3-6J^{OveA1#;5cRlIfK^ME_}Z-F-w0T(uc@1MbjJkYQ~cVV&}lvm zLa#e9ejSR@)_rb?mv|9mm9RMQ1)c$684#NRp#@$NKf-5dn=p-DniK=Zl6Vp6r5)m3 zq+9p=4OXo=tV4W>oX$nz);+o+zgy2|KYZZn^;*TAU#{M-$R1@ToM70Po)5!fIK^)( zHS?hu6$+gb-`YQmEi)y4n}1plb%FRB>(;Vx(vFLh9b1}4<5b%bW_&7`$9T$;Vj+%Y oJ&oTU#6J$1mtjeVs3jq`4a2BIohhkwk?cGdQG?B&U(t-6e}`E4NB{r; diff --git a/controller.py b/controller.py index 1e2c2a2..312d239 100644 --- a/controller.py +++ b/controller.py @@ -66,6 +66,10 @@ def settlement_prediction(business_code, cons_code): if row[0] is None or row[1] is None or row[2] is None: continue + # [로그 추가] 침하량이 0인 데이터 찾기 + if float(row[0]) == 0: + print(f"[DB_CHECK] Zero settlement found - Business: {business_code}, Cons: {cons_code}, Day(nod): {row[2]}") + settlement.append(float(row[0])) surcharge.append(float(row[1])) time.append(float(row[2])) diff --git a/settle_prediction_steps_main.py b/settle_prediction_steps_main.py index b051f45..83bc389 100644 --- a/settle_prediction_steps_main.py +++ b/settle_prediction_steps_main.py @@ -25,17 +25,29 @@ from scipy.interpolate import interp1d # Function 섹션 # ================= -# 주어진 계수를 이용하여 쌍곡선 시간-침하 곡선 반환 +# 주어진 계수를 이용하여 쌍곡선 시간-침하 곡선 반환 # thkim 0으로 나누기 오류 수정 def generate_data_hyper(px, pt): - return pt / (px[0] * pt + px[1]) + denom = px[0] * pt + px[1] + denom_safe = np.where(denom == 0, 1e-10, denom) + return pt / denom_safe -# 주어진 계수를 이용하여 아사오카 시간-침하 곡선 반환 +# 주어진 계수를 이용하여 아사오카 시간-침하 곡선 반환 # thkim 0으로 나누기 오류 수정 def generate_data_asaoka(px, pt, dt): - return (px[1] / (1 - px[0])) * (1 - (px[0] ** (pt / dt))) + # (1 - px[0])이 0이 되는 경우 방지 + denom = 1 - px[0] + if denom == 0: + denom = 1e-10 + + # dt(간격)가 0인 경우 방지 + dt_safe = dt if dt != 0 else 1e-10 + + return (px[1] / denom) * (1 - (px[0] ** (pt / dt_safe))) # 회귀식과 측정치와의 잔차 반환 (비선형 쌍곡선) -def fun_hyper_nonlinear(px, pt, py): - return pt / (px[0] * pt + px[1]) - py +def fun_hyper_nonlinear(px, pt, py): # thkim 0으로 나누기 오류 수정 + denom = px[0] * pt + px[1] + denom_safe = np.where(denom == 0, 1e-10, denom) + return pt / denom_safe - py # 회귀식과 측정치와의 잔차 반환 (가중 비선형 쌍곡선) def fun_hyper_weight_nonlinear(px, pt, py, pw): @@ -43,6 +55,10 @@ def fun_hyper_weight_nonlinear(px, pt, py, pw): # 회귀식과 측정치와의 잔차 반환 (기존 쌍곡선) def fun_hyper_original(px, pt, py): + # 나누기 에러 발생 시 로그 출력 + if np.any(py == 0): + print(f"[DEBUG] py(settlement) contains zero value! pt sample: {pt[:5]}") + py_safe = np.where(py == 0, 1e-10, py) # thkim 0으로 나누기 오류 수정 return px[0] * pt + px[1] - pt / py # 회귀식과 측정치와의 잔차 반환 (아사오카) @@ -298,8 +314,33 @@ def run_settle_prediction(point_name, np_time, np_surcharge, np_settlement, # 회귀분석 시행 (기존 쌍곡선법) - (0, 0)에 해당하는 초기 데이터를 제외하고 회귀분석 실시 x0 = np.ones(2) - res_lsq_hyper_original = least_squares(fun_hyper_original, x0, - args=(tm_hyper[1:], sm_hyper[1:])) + # [로그 추가] 입력 데이터의 상태를 확인 + print(f"[LOG] {point_name} - tm_hyper size: {len(tm_hyper)}, sm_hyper contains zero: {np.any(sm_hyper == 0)}") + if np.any(sm_hyper <= 0): + print(f"[LOG] Problematic sm_hyper values: {sm_hyper[sm_hyper <= 0]}") + + # 1. 침하량이 0보다 큰 유효한 데이터만 필터링 (0으로 나누기 근본적 방지) + valid_indices = np.where(sm_hyper > 0)[0] + + # 2. 필터링된 데이터가 최소 2개 이상인지 확인 (회귀분석을 위한 최소 조건) + if len(valid_indices) >= 2: + tm_hyper_valid = tm_hyper[valid_indices] + sm_hyper_valid = sm_hyper[valid_indices] + + # (0, 0) 제외 로직이 이미 필터링에 포함되어 있으므로 [1:] 없이 수행 가능 + x0 = np.ones(2) + try: + res_lsq_hyper_original = least_squares(fun_hyper_original, x0, + args=(tm_hyper_valid, sm_hyper_valid)) + x_hyper_original = res_lsq_hyper_original.x + except ValueError as e: + print(f"[Warning] Optimization failed for Original Hyperbolic: {e}") + x_hyper_original = np.ones(2) * 0.001 # 실패 시 기본값 세팅 + else: + # 데이터가 부족할 경우 에러 대신 기본값이나 알림 처리 + print(f"[Warning] {point_name}: Not enough valid settlement data (>0) for Original Hyperbolic.") + x_hyper_original = np.ones(2) * 0.001 + # 기존 쌍곡선 법 계수 저장 및 출력 x_hyper_original = res_lsq_hyper_original.x