S660に続き、アルトワークスのCANを解析して、馬力を計測してみた。S660の場合は、HondaのCAN解析を行っているプロジェクトがあったのでそちらを参考にしたが、HA36Sアルトワークスの場合はCANデータ解析の情報がなかったので、自分で解析を行う必要があった。車速データ探しから始め、馬力計算できるようになったので記事にしてみた。
S660の馬力計測の記事はこちら
Arduino＋CAN-BUS Shieldを用いてCAN受信プログラムを作成し、SDカードに保存。保存した速度データを後から解析することで馬力を計算している。今回の手法を用いれば、何か車両に変更があった時にそれがパワーにどういう影響があるのかを考察することができる。
馬力計算の原理
馬力＝トルク×回転数で算出が可能なことを利用し、車両からCAN経由で取得した速度の時系列データから加速度を計算。加速度からトルクと馬力を算出する。
馬力とトルクの関係はこちらの記事を参照
トルクの算出
詳細の物理については知見が深いサイトに説明を譲るとして、トルクは、
で算出することができる。
トルクを算出する情報のうち、変速比・ファイナルギア比・車重はMT車ならカタログに載っている値で計算できる。アルトワークスの場合は、オートマの場合でもCVTではなく5AGSが採用されているため、この手法はMT/AT問わず使用可能。
タイヤ半径については、アルトワークスの場合は前後同径。今回は駆動輪が前輪なので、前輪のサイズを用いることとした。加速度は時間あたりの速度増加量なので、速度の時系列データを収集した。
馬力を計算する為に必要となる回転数は、ギア固定（2速）で走行。
速度から回転数を算出する手法についてはこちらの記事を参照
ということで、速度の時系列データを取得して、トルクを算出。速度から回転数を計算し、トルク×回転数で馬力を算出してみた。
速度データの取得（Arduino＋CAN-BUS Shield V2.0）
ArduinoでCANデータを取得するモジュールは作成済み。今回はこちらを用いてCANデータを取得する。
Arduinoで取得するCANデータIDの選定
スズキのアルトワークスのCAN解析結果は、Webページを探しても出てこなかった。Twitterで0x314に速度データが流れているとの情報があり、これをきっかけに、アルトワークスのCAN解析を行った。0x1B8にESP(エレクトロニック・スタビリティ・プログラム)と思われる車輪速のデータが流れていることを突き止めた。馬力計算のデータは、周期が速い方が確度が上がる。今回は20ms周期で出力されていると思われる0x1B8を使うこととした。
0x1B8に流れているESP車輪速情報
それぞれビッグエンディアン、符号なし。
Byte0,Byte1のように2バイトを使用して16バイトで一つの車輪速を表している。係数は0.025。単位はkm/h。
Byte0,1は右前輪、Byte2,3は左前輪、Byte4,5は右後輪、Byte6,7は左後輪。
これらのうち、前輪の左右輪情報を使って馬力を計測してみた。フィルターでID 0x1B8のみ取得してSDカードに保存（Arduinoプログラムについてはこの記事の後半に記載）。
実走でのデータ取得
最大出力を計算するには、平坦な場所でアクセル全開のフル加速時のデータを取得する必要がある。HA36S(MT)の2速で7,000回転時の速度は75km/hとなるため、フル加速データを安全に計測できる場所は自ずと限られてくる。実際に取得を考えている場合は法的にも、安全にも配慮して計測して欲しい。
ArduinoとCAN-BUS Shield、OBD2接続コネクタの作成
CANを車体から取得するのはOBD2カプラが便利。こんな感じで、S660のときに作成したArduino＋CAN-BUS ShieldにOBD2コネクタを接続した。Arduinoの電源はスマホ充電用のUSBを使用。
アルトワークスのOBD2コネクタの場所
車体側のOBD2コネクタは運転席のアクセルペダルの上の方にある。
今回はここに自作のCANケーブルを接続し、ArduinoのCAN-BUS Shieldに接続した。
ArduinoのSDカードに保存されるデータフォーマット
取得したデータは、時系列のCSV形式とした。Arduino側で、以下フォーマットでSDカードに格納されるようにプログラムを作成。Timestamp[ms], CAN ID(dec), DataLength, Byte0(dec), Byte1(dec), …. Byte7(dec) でSDカードに保存される。実際にArduinoで取得したデータの抜粋を下記に示す。
19680,440,8,3,60,3,59,3,49,3,47
19692,440,8,3,63,3,63,3,52,3,51
19704,440,8,3,66,3,65,3,57,3,53
19716,440,8,3,71,3,69,3,61,3,57
19728,440,8,3,74,3,72,3,63,3,60
19741,440,8,3,78,3,77,3,67,3,64
19753,440,8,3,81,3,79,3,69,3,67
19765,440,8,3,84,3,83,3,74,3,70
19777,440,8,3,89,3,87,3,77,3,74
19789,440,8,3,93,3,90,3,82,3,77
19801,440,8,3,95,3,96,3,82,3,81
19813,440,8,3,97,3,101,3,86,3,85
19825,440,8,3,104,3,103,3,89,3,87
19837,440,8,3,109,3,107,3,94,3,92CANデータの解析と、グラフ表示
SDカードで取得したデータを解析。CANデータの解析方法については、こちらの記事で解説。
解析した結果、今回取得したデータをグラフにしたのがこちら。
7回の全力加速（2速）データが取得できている。この速度データから
・サンプリング時間あたりの速度差を計算し、加速度を算出する。
・速度に対応した回転数を算出する。
そして、加速度・車重・ギア比・タイヤ半径からトルクを計算し、対応する回転数から馬力を計算するプログラムを作成した（解析用ソースコードはこの記事の最後に添付）。
HA36Sのパワーカーブ
取得したCAN速度データをもとに、Pythonで解析を実施。以下グラフは、横軸回転数、縦軸はトルクと馬力を示している。
今回の計測では、最大トルクは3,500rpmで78N･m。最大出力は4,700rpmで46馬力程度と出ている。エンジン出力で記載されているカタログの値（最大トルク発生は3,000rpm、最大出力は6,000rpm）とは異なる理由は、空気抵抗やタイヤの転がり抵抗、慣性抵抗を含んだ馬力となっているためだろう。特に、空気抵抗は速度の2乗に比例して大きくなるため、高回転領域では、（見かけ上の）出力が低下していると考えられる。
3,000rpm付近は慣性抵抗、6,000rpm付近は空気抵抗による影響を受けているのだろうと思われる。空気抵抗や転がり抵抗は計算によって算出できるので補正は可能だが、今回は行っていない。自分の車の特性把握や、パーツを交換したときの効果を測定するには、メーカーが謳うエンジン単体の性能ではなく「リアルワールド」での特性を把握したかったからだ。
アルトワークスのエンジンの特徴
今回算出したグラフから、アルトワークスのトルクは低回転から出力を大きくする典型的な小排気量ターボ車であることがわかる。S660の場合もそうであったが5,000rpm以降はトルクが出にくい特性だ。低回転での「キビキビ」さを重視したセッティングとなっている事がわかる。また、後半の馬力は頭打ち。スロットルをわざと出力を絞っているような形跡がある（おそらく64馬力規制のためだろう）。
今回のCAN解析を行う中で、0x13Fにエンジン回転数とスロットル開度(？)が出力されていることを発見したので、そちらのデータも併せてみてみることにした。
下のグラフは全開加速時のデータ。↔の区間は、アクセルペダルを床まで踏み込んでいるが、スロットル開度（もしくは要求トルク、あるいはブースト？？）と思われる信号は、アクセル全開にも関わらず回転数の上昇にともなって減り続けている（赤い◯の部分）。
また、初期応答の反応を鈍らせるような「なまし」の制御が入っており（黄色◯部分）、低速時に扱いやすい特性である一方、レスポンスの悪さを感じる。これは、急なトルク変化でホイールスピンを避けるためのセッティングによるものだろう。Suzukiのエンジニアが思慮を重ねた結果の設計のため、当然様々なノウハウがあるのだろうが、一部のユーザーにはもう少しハイレスポンスのマップがあっても良さそうだ（ただしその分、街乗りで乗りにくくなるだろうが）。
このグラフと、パワーカーブを比べてみた結果が下の図。スロットル制御と思われる制御波形と、トルク波形が同じような形状であることから、アルトワークスの高回転領域のトルク減少は、回転数に応じてスロットルを絞り込む事による制御の影響を受けている事が想像できる。おそらく、64馬力の規制に対応するための制御であると考えられる。個人的には、軽自動車の自主規制のために、色々苦労してるんだな～という感想。
以上を踏まえて、エンジンチューニングを行うなら、まずは4,000rpm以降で（おそらく）わざと絞っているスロットルを開けるような制御を入れるのが常道かと思われる。ただ、最近の車両がどこまでの情報をもとにエンジン制御しているのか不明なので、ブローが怖い。この辺は経験豊富なチューニングショップであれば、Suzukiのエンジニアが意図的に絞り込んだ領域を、マージンを確保しながら開放する事ができるはずだ。
一方で、アクセルレスポンス（いわゆる「ツキ」）は、もう少し敏感なセッティングにした方が私好みだ。アクセル開度に対するスロットル開度の追従性をプログラム的に上げることができればかなり改善すると思われる。
今回使用したデバイス
今回、データ取得のために使用したデバイスはこちら。その他銅線を含めても、1万円程度で揃う。Arduino Unoを使ったり、互換性のある中華製Arduinoを使ったり等すればもっと安くできるはずだ。今回はS660の解析の際に使っていたものをそのまま流用。
Arduino
Arduino Unoでも良いと思うが、色々遊んでみたかったので、Arduino megaを選定。
CAN-BUS Shield V2.0
OBDコネクタ
D-SUB9ピン メスコネクタ
CAN-BUS Shield V2.0と車両を接続する際に必要なコネクタ。はんだ付けして使う。
データ取得プログラム（Arduino）
ArduinoとCAN-BUS Shield を組み合わせて車両からCANデータを取得し、実装されているmicroSDカードスロットルを利用して、microSDカードに取得したCANデータを保存するプログラムを作成した（こちら）。
データ解析プログラム（Python）
岡本一車プロのレビュアーや解析のエキスパートの皆様、算出の考え方や計算方法等、アドバイスがありましたらコメントお願いします。
"""
File: power_analysis_HA36S.py
Author: 岡本一車
https://kurumashikou.com/
"""
#必要なライブラリのimport
import numpy as np;
import csv
import math;
from matplotlib import pyplot as plt;
# ログファイルの定義 ----------------------------------------------------------------------------------------
# プログラムディレクトリへのパス
dirPath = "D:power_analysis";
# ログディレクトリへのパス
logFileDir = dirPath + "\\ha36s_candata";
# ログファイルへのパス
logFilePath = logFileDir + "\\" + "DATALOG_.CSV";
#車両
car_name = "HA36S";
""" 車両個別設定項目 ----------------------------------------------------------------------------------------
車重
タイヤ
ギア比
を設定する。
------------------------------------------------------------------------------------------------------------ """
# 車重[kg]
weight = 670 + 100;
# タイヤの設定 -----------------------------------------
# タイヤ太さ[mm] ex:165/55R15なら165
tire_width_mm = 165;
# 扁平率[%] ex:165/55R15なら55
tire_aspect_ratio = 55;
# ホイールサイズ[inch] ex:165/55R15なら15
wheel_size = 15;
# ギアレシオ -------------------------------------------
"""
HA36Sの場合
final: 4.705
1st: 3.545
2nd: 2.105
3rd: 1.521
4th: 1.148
5th: 0.897
"""
gear_ratio_arr = [
4.705,
3.545,
2.105,
1.521,
1.148,
0.897
];
# 計測に使用したギア
gear_select = 2;
# 走行するギアレシオを設定する
gear_ratio = gear_ratio_arr[gear_select];
final_gear_ratio = gear_ratio_arr[0];
# データ取得の最低回転数
min_rpm = 2500;
max_rpm = 7000;
""" ----------------------------------------------------------------------------------------------------------
移動平均を算出する関数
arg1:データのlist
num:移動平均算出の個数
return 移動平均算出後のLIST
データ個数が不足している場合、0が入る
例えば… moving_Average([2, 4, 9, 17],3) の結果、
returnは[0, 0, 5.0, 10.0]
(2+4+9)/3 ↑ ↑(4+9+17)/3
---------------------------------------------------------------------------------------------------------- """
def moving_Average(in_list, num):
cumsum = [0];
ret_mvingAve = [];
for i, x in enumerate(in_list, 1):
cumsum.append(cumsum[i-1] + x);
if i>=num:
elem_movingAve = (cumsum[i] - cumsum[i-num])/num;
ret_mvingAve.append(elem_movingAve);
for j in range(num-1):
ret_mvingAve.insert(0,0);
return ret_mvingAve;
""" ----------------------------------------------------------------------------------------------------------
ログデータの読み込み処理
読み込むCSVファイルフォーマット
TimeStamp[ms],CAN-ID,DLC, ...
DateByte1, DateByte2, DateByte3, DateByte4, DateByte5, DateByte6, DateByte7, DateByte8
---------------------------------------------------------------------------------------------------------- """
# CSV形式のログファイルから、データを読み込む
csv_file = open( logFilePath, "r", encoding="ms932", errors="", newline="");
# CSVファイルからデータを読み込むためのReaderオブジェクトを設定
CanLogReader = csv.reader(csv_file);
# 生のログ全データを抽出し、保持
CANLogRawData = [row for row in CanLogReader];
# NumPyのArray型に変更
CANLogRawData = np.array(CANLogRawData);
# CSVファイルのClose
csv_file.close();
""" ----------------------------------------------------------------------------------------------------------
Listデータの型を数値型に変換する
-------------------------------------------------------------------------------------------------------- """
for i in range(len(CANLogRawData)):
for j in range(len(CANLogRawData[i])):
# 格納されている値が10進数かどうかを確認して、10進数に変換できないデータは-1に変換する
if(not CANLogRawData[i][j].isdecimal()):
CANLogRawData[i][j] = -1;
# 全体をキャストしてint型に変更
CANLogRawData = CANLogRawData.astype(np.int_);
""" ----------------------------------------------------------------------------------------------------------
結果算出用に、時間およびΔtを算出する
-------------------------------------------------------------------------------------------------------- """
# 時間の取得
log_time = np.array(CANLogRawData[:,0]);
# 差分取得用に、1サンプル進んだt1を作成
delta_t1 = np.array(log_time[1:]);
# log_time
delta_t = delta_t1 - log_time[0:-1];
# 結果の先頭に0を追加
delta_t = np.append([0], delta_t);
# CANデータ部だけを抽出する
canData = np.array(CANLogRawData[:,3:11]);
#FR車速の取得 --------------------------------------------------------------------------------------------
# |byte0 |byte1 |byte2 |byte3 |byte4 |byte5 |byte6 |byte7 |
# |11111111|11111111|00000000|00000000|00000000|00000000|00000000|00000000|
# |<-----data ----->|
# --------------------------------------------------------------------------------------------------------
# 先頭からbyte0,byte1を抜き出す
WHEEL_SPEED_FR = canData[:,0:2];
# 抜き出したbyteに対して、対象のビットのみ取得する演算
WHEEL_SPEED_FR = (WHEEL_SPEED_FR & [0b11111111,0b11111111]);
# byte1の最下位ビットが１桁目となるようにビットシフト、byte0の最下位ビットが8桁目となるようにシフト
WHEEL_SPEED_FR = ([WHEEL_SPEED_FR[:,0] << 8 , WHEEL_SPEED_FR[:,1]]);
# byte0とbyte1を足し合わせて、Factor0.025をかける
WHEEL_SPEED_FR = np.sum(WHEEL_SPEED_FR, axis=0)*0.025;
#FL車速の取得 --------------------------------------------------------------------------------------------
# |byte0 |byte1 |byte2 |byte3 |byte4 |byte5 |byte6 |byte7 |
# |00000000|00000000|11111111|11111111|00000000|00000000|00000000|00000000|
# | |<-----data ----->|
# --------------------------------------------------------------------------------------------------------
WHEEL_SPEED_FL = canData[:,2:4];
WHEEL_SPEED_FL = (WHEEL_SPEED_FL & [0b11111111,0b11111111]);
WHEEL_SPEED_FL = ([WHEEL_SPEED_FL[:,0] << 8 , WHEEL_SPEED_FL[:,1]]);
WHEEL_SPEED_FL = np.sum(WHEEL_SPEED_FL, axis=0)*0.025;
#RR車速の取得 --------------------------------------------------------------------------------------------
# |byte0 |byte1 |byte2 |byte3 |byte4 |byte5 |byte6 |byte7 |
# |00000000|00000000|00000000|00000000|11111111|11111111|00000000|00000000|
# |<-----data ----->|
# --------------------------------------------------------------------------------------------------------
WHEEL_SPEED_RR = canData[:,4:6];
WHEEL_SPEED_RR = (WHEEL_SPEED_RR & [0b11111111,0b11111111]);
WHEEL_SPEED_RR = ([WHEEL_SPEED_RR[:,0] << 8 , WHEEL_SPEED_RR[:,1]]);
WHEEL_SPEED_RR = np.sum(WHEEL_SPEED_RR, axis=0)*0.025;
#RL車速の取得 --------------------------------------------------------------------------------------------
# |byte0 |byte1 |byte2 |byte3 |byte4 |byte5 |byte6 |byte7 |
# |00000000|00000000|00000000|00000000|00000000|00000000|11111111|11111111|
# | |<-----data ----->|
# --------------------------------------------------------------------------------------------------------
WHEEL_SPEED_RL = canData[:,6:8];
WHEEL_SPEED_RL = (WHEEL_SPEED_RL & [0b11111111,0b11111111]);
WHEEL_SPEED_RL = ([WHEEL_SPEED_RL[:,0] << 8 , WHEEL_SPEED_RL[:,1]]);
WHEEL_SPEED_RL = np.sum(WHEEL_SPEED_RL, axis=0)*0.025;
# リア車速の平均を求め、馬力算出用の速度とする
vel_kph = np.array([[WHEEL_SPEED_FL],[WHEEL_SPEED_FR]]);
vel_kph = vel_kph.mean(axis=(0,1));
# 算出された速度に移動平均フィルタをかける
vel_kph = np.array(moving_Average(vel_kph,4));
"""----------------------------------------------------------------------------------------------------------
速度[km/h]から加速度[m/s2]を算出
----------------------------------------------------------------------------------------------------------"""
# 速度の差分を計算
delta_v1 = np.array(vel_kph[1:]);
# 1データ前の速度との差分を取得
delta_v = delta_v1 - vel_kph[0:-1];
# 1データあたりの加速度を算出し、サンプリングタイムからm/s2に変換
accel_mps2 = delta_v/3.6*1000/delta_t[1:];
# 結果の先頭に0を追加（1番目のデータは0とする）
accel_mps2 = np.append([0], accel_mps2);
# 加速度に移動平均フィルタをかける。加速度は重めのフィルタが必要
accel_mps2 = np.array(moving_Average(accel_mps2,15));
"""----------------------------------------------------------------------------------------------------------
速度[km/h]から回転数[rpm]を算出
----------------------------------------------------------------------------------------------------------"""
# タイヤ直径[m] = ((2 × タイヤ太さ[mm] × 扁平率[%]/100) + (ホイールサイズ[inch] × 25.4(インチ変換定数))/1000
tire_diameter = ((2*tire_width_mm*tire_aspect_ratio/100)+(wheel_size*25.4))/1000
# タイヤ外周[m] = タイヤ直径 × π
tire_circumference = tire_diameter*math.pi;
# 回転数[rpm] = (速度[km/h]×ギア比×ファイナルギア比)/タイヤ外周[m]×1000/60
tacho_rpm = np.array((vel_kph*gear_ratio*final_gear_ratio)/tire_circumference*1000/60);
"""----------------------------------------------------------------------------------------------------------
加速度[m/s2]からトルク[N･m]を算出
----------------------------------------------------------------------------------------------------------"""
# トルク[N･m] = 加速度[m/s2]/(ギア比×ファイナル) × タイヤ直径/2 × 車重[kg]
torq_nm = np.array(accel_mps2/(gear_ratio*final_gear_ratio)*(tire_diameter/2)*weight);
"""----------------------------------------------------------------------------------------------------------
トルク[N･m]と回転数からパワー[kw]と馬力[ps]を算出
----------------------------------------------------------------------------------------------------------"""
# パワー[Kw] = 2π x トルク(N・m) x 回転数(rpm) / 60 / 1000
power_kw = np.array(2*math.pi*torq_nm*tacho_rpm/60/1000);
# トルク、パワーの計算結果を格納(回転数、トルク、パワー、馬力)
power_data = np.array([tacho_rpm, torq_nm, power_kw, power_kw*1.35962]);
"""----------------------------------------------------------------------------------------------------------
取得したデータを整形する。不要なデータの削除およびマイナス値の削除、並び替え
----------------------------------------------------------------------------------------------------------"""
# 数値マイナスのものを削除
power_data = np.delete(power_data, np.where(power_data<0)[1], axis=1);
# 回転数がmin_rpm～max_rpmで指定した範囲に絞り込み
power_data = power_data[:,(power_data[0]>min_rpm)&(power_data[0]<max_rpm)];
# 馬力データを回転数で並び替え
sorted_index = np.argsort(power_data[0]);
power_data = power_data[:,sorted_index];
# ここまでの処理で、以下の配列が得られた
# power_data =
# [
# [昇順にソートされた回転数rpmの配列],
# [回転数に対応したトルクN･m],
# [回転数に対応したパワーkw],
# [回転数に対応した馬力ps]
# ]
# ----------------------------------------------------------------------------------------------------------
"""----------------------------------------------------------------------------------------------------------
近似曲線の作成
----------------------------------------------------------------------------------------------------------"""
# 近似曲線の始点、終点の定義。最高回転数のレブリミッターが効くので、そのより少し低めに設定
poly_x_ax = list(range(min_rpm,max_rpm-200,10));
poly_y_axTorq = np.polyfit(power_data[0,:], power_data[1,:], 10);
poly_y_axPS = np.polyfit(power_data[0,:], power_data[3,:], 10);
"""----------------------------------------------------------------------------------------------------------
グラフの描写
----------------------------------------------------------------------------------------------------------"""
# Figureの作成
velfig = plt.figure();
velfig_ax1 = velfig.add_subplot(111);
velfig_ax1.set_title("vel kph " + car_name);
velfig_ax1.set_xlabel("time[ms]");
velfig_ax1.set_ylabel("vel[kph]");
plt.plot(log_time,vel_kph);
# Figureの作成
fig = plt.figure();
# 軸を追加
ax1 = fig.add_subplot(111);
# 軸の設定（回転数、トルク）
ax1.set_xlim(min_rpm-100, max_rpm+100);
ax1.set_ylim(0, int(max(power_data[1,:])*1.2));
ax1.set_title("Power and Torque curve " + car_name);
ax1.set_xlabel("EnginSpeed[rpm]");
ax1.set_ylabel("torque[N.m]");
# トルクの散布図を追加
plt.scatter(power_data[0,:], power_data[1,:], s=1, color="cornflowerblue");
# トルクの近似曲線の追加
ax1.plot(poly_x_ax, np.poly1d(poly_y_axTorq)(poly_x_ax), color="blue", label="Torque");
# 軸の設定（馬力側）
ax2 = ax1.twinx();
ax2.set_ylim(0, int(max(power_data[3,:])*1.2));
ax2.set_ylabel("Power[PS]");
# 馬力の散布図を追加
ax2.scatter(power_data[0,:], power_data[3,:], s=1, color="lightsalmon");
# 馬力の近似曲線の追加
ax2.plot(poly_x_ax, np.poly1d(poly_y_axPS)(poly_x_ax), color="red", label="Power[PS]");
# グリッド線の追加
ax1.grid(b=None, which="both", axis="x");
ax2.grid(b=None, which="both", axis="both");
# ラベルの追加
h1, l1 = ax1.get_legend_handles_labels();
h2, l2 = ax2.get_legend_handles_labels();
ax2.legend(h1+h2, l1+l2, loc="lower right");
# グラフを描画
plt.show();
コメント