##############################################################
## Calculate Vibration Transmission (dB)
##############################################################

jobName = "Name of Job"
stepName = "Name of Steady State Dynamic Step"
outputSetName = "Name of Output Point Set"

from odbAccess import*
from abaqusConstants import*
import string
import numpy as np
import os

odb = openOdb(path = jobName+'.odb')

outfile = open(jobName + '.csv', 'w')
outfile.write('Frequency [Hz]' + ',' + 'Transmission [dB]' + '\n')

for fm in range(1, len(odb.steps[stepName].frames)):

  timeFrame = odb.steps[stepName].frames[fm]
  readNode = odb.rootAssembly.nodeSets[outputSetName.upper()]
  Disp = timeFrame.fieldOutputs['U']
  readNodeDisp = Disp.getSubset(region=readNode)
  readNodeDispValues = readNodeDisp.values

  frequency = np.zeros(len(odb.steps[stepName].frames))
  Transmission = np.zeros(len(odb.steps[stepName].frames))
  Disp_real = np.zeros(len(odb.steps[stepName].frames))
  Disp_imaginary = np.zeros(len(odb.steps[stepName].frames))
  Disp_magnitude = np.zeros(len(odb.steps[stepName].frames))
  frequency[fm] = timeFrame.frameValue
  Disp_real[fm] = readNodeDispValues[0].data[2] # 0-X Direction; 1-Y Direction; 2-Z Direction
  Disp_imaginary[fm] = readNodeDispValues[0].conjugateData[2]
  Disp_magnitude[fm] = sqrt(Disp_real[fm]**2.0 +Disp_imaginary[fm]**2.0)
  Transmission[fm] = 20.0*log10(Disp_magnitude[fm])

  outfile.write(str(frequency[fm]) + ',' + str(Transmission[fm]) + ',' + '\n')

outfile.close()

odb.close()