211 lines
3.8 KiB
Go
211 lines
3.8 KiB
Go
package analize
|
|
|
|
import (
|
|
"encoding/csv"
|
|
"errors"
|
|
"fmt"
|
|
"math"
|
|
"os"
|
|
"strconv"
|
|
"strings"
|
|
|
|
"gonum.org/v1/gonum/dsp/fourier"
|
|
"gorm.io/gorm"
|
|
)
|
|
|
|
func AnalyzeCSV(db *gorm.DB, caseID uint, caseName, csvPath string) (*Analize, error) {
|
|
file, err := os.Open(csvPath)
|
|
if err != nil {
|
|
return nil, err
|
|
}
|
|
defer file.Close()
|
|
|
|
reader := csv.NewReader(file)
|
|
records, err := reader.ReadAll()
|
|
if err != nil {
|
|
return nil, err
|
|
}
|
|
if len(records) < 2 {
|
|
return nil, fmt.Errorf("csv has no data rows")
|
|
}
|
|
|
|
headers := make(map[string]int, len(records[0]))
|
|
for i, header := range records[0] {
|
|
headers[strings.ToLower(strings.TrimSpace(header))] = i
|
|
}
|
|
|
|
getIndex := func(name string) (int, error) {
|
|
idx, ok := headers[strings.ToLower(name)]
|
|
if !ok {
|
|
return -1, fmt.Errorf("missing column %q", name)
|
|
}
|
|
return idx, nil
|
|
}
|
|
|
|
tIndex, err := getIndex("t")
|
|
if err != nil {
|
|
return nil, err
|
|
}
|
|
psiMIndex, err := getIndex("psi_m")
|
|
if err != nil {
|
|
return nil, err
|
|
}
|
|
psiLIndex, err := getIndex("psi_l")
|
|
if err != nil {
|
|
return nil, err
|
|
}
|
|
|
|
times := make([]float64, 0, len(records)-1)
|
|
psiM := make([]float64, 0, len(records)-1)
|
|
psiL := make([]float64, 0, len(records)-1)
|
|
|
|
for _, record := range records[1:] {
|
|
if len(record) <= maxInt(tIndex, psiMIndex, psiLIndex) {
|
|
continue
|
|
}
|
|
|
|
t, err := parseFloat(record[tIndex])
|
|
if err != nil {
|
|
continue
|
|
}
|
|
m, err := parseFloat(record[psiMIndex])
|
|
if err != nil {
|
|
continue
|
|
}
|
|
l, err := parseFloat(record[psiLIndex])
|
|
if err != nil {
|
|
continue
|
|
}
|
|
|
|
times = append(times, t)
|
|
psiM = append(psiM, m)
|
|
psiL = append(psiL, l)
|
|
}
|
|
|
|
if len(psiL) == 0 {
|
|
return nil, fmt.Errorf("no numeric rows found in csv")
|
|
}
|
|
|
|
psiMax := maxFloat(psiM)
|
|
psiLMax := maxFloat(psiL)
|
|
omega := mainAngularFrequency(times, psiL)
|
|
|
|
result := &Analize{
|
|
Name: fmt.Sprintf("analysis-%d", caseID),
|
|
CaseID: caseID,
|
|
CaseName: caseName,
|
|
PsiMax: psiMax,
|
|
PsiLMax: psiLMax,
|
|
Omega: omega,
|
|
}
|
|
|
|
var existing Analize
|
|
err = db.Where("case_id = ?", caseID).First(&existing).Error
|
|
if err != nil {
|
|
if errors.Is(err, gorm.ErrRecordNotFound) {
|
|
if err := db.Create(result).Error; err != nil {
|
|
return nil, err
|
|
}
|
|
return result, nil
|
|
}
|
|
return nil, err
|
|
}
|
|
|
|
result.ID = existing.ID
|
|
if err := db.Save(result).Error; err != nil {
|
|
return nil, err
|
|
}
|
|
return result, nil
|
|
}
|
|
|
|
func parseFloat(value string) (float64, error) {
|
|
return strconv.ParseFloat(strings.TrimSpace(value), 64)
|
|
}
|
|
|
|
func maxInt(values ...int) int {
|
|
max := 0
|
|
for _, v := range values {
|
|
if v > max {
|
|
max = v
|
|
}
|
|
}
|
|
return max
|
|
}
|
|
|
|
func maxFloat(values []float64) float64 {
|
|
if len(values) == 0 {
|
|
return 0
|
|
}
|
|
max := values[0]
|
|
for _, v := range values[1:] {
|
|
if v > max {
|
|
max = v
|
|
}
|
|
}
|
|
return max
|
|
}
|
|
|
|
func mainAngularFrequency(times, values []float64) float64 {
|
|
if len(values) < 2 || len(times) < 2 {
|
|
return 0
|
|
}
|
|
|
|
dt := averageDelta(times)
|
|
if dt <= 0 {
|
|
return 0
|
|
}
|
|
|
|
centered := make([]float64, len(values))
|
|
mean := 0.0
|
|
for _, v := range values {
|
|
mean += v
|
|
}
|
|
mean /= float64(len(values))
|
|
for i, v := range values {
|
|
centered[i] = v - mean
|
|
}
|
|
|
|
fft := fourier.NewFFT(len(centered))
|
|
coeffs := fft.Coefficients(nil, centered)
|
|
if len(coeffs) < 2 {
|
|
return 0
|
|
}
|
|
|
|
bestIndex := 1
|
|
bestAmp := 0.0
|
|
for i := 1; i < len(coeffs); i++ {
|
|
amp := absComplex(coeffs[i])
|
|
if amp > bestAmp {
|
|
bestAmp = amp
|
|
bestIndex = i
|
|
}
|
|
}
|
|
|
|
frequencyHz := fft.Freq(bestIndex) / dt
|
|
return 2 * math.Pi * frequencyHz
|
|
}
|
|
|
|
func absComplex(v complex128) float64 {
|
|
return math.Hypot(real(v), imag(v))
|
|
}
|
|
|
|
func averageDelta(times []float64) float64 {
|
|
if len(times) < 2 {
|
|
return 0
|
|
}
|
|
total := 0.0
|
|
count := 0
|
|
for i := 1; i < len(times); i++ {
|
|
delta := times[i] - times[i-1]
|
|
if delta <= 0 {
|
|
continue
|
|
}
|
|
total += delta
|
|
count++
|
|
}
|
|
if count == 0 {
|
|
return 0
|
|
}
|
|
return total / float64(count)
|
|
}
|