In [22]:
import tensorflow as tf
import numpy
import scipy.io
import pandas
import os
from tensorflow.python.client import timeline
import time
#carico file dati
percorsoFile = "/home/protoss/Documenti/TESI/DATI/peaks.mat"
#print(picchi.shape)
#picchi[0]
#nb: picchi ha 0-tempi
# 1-frequenze
# 4-pesi
#ora popolo il dataframe
tabella = pandas.DataFrame(scipy.io.loadmat(percorsoFile)['PEAKS'])
tabella.drop(tabella.columns[[2, 3]], axis = 1, inplace=True)
tabella.columns = ["tempi", "frequenze","pesi"]
#fascia di sicurezza
securbelt = 4000
headerFreq= scipy.io.loadmat(percorsoFile)['hm_job'][0,0]['fr'][0]
headerSpindown = scipy.io.loadmat(percorsoFile)['hm_job'][0,0]['sd'][0]
epoca = scipy.io.loadmat(percorsoFile)['basic_info'][0,0]['epoch'][0,0]
#nb: headerFreq ha 0- freq minima,
# 1- step frequenza,
# 2- enhancement in risoluzone freq,
# 3- freq massima,
#headerSpindown ha 0- spin down iniziale di pulsar
# 1- step spindown
# 2- numero di step di spindown
#Definisco relative variabili per comodità e chiarezza del codice
#frequenze
minFreq = headerFreq[0]
maxFreq = headerFreq[3]
enhancement = headerFreq[2]
stepFrequenza = headerFreq[1]
stepFreqRaffinato = stepFrequenza/enhancement
freqIniz = minFreq- stepFrequenza/2 - stepFreqRaffinato
freqFin = maxFreq + stepFrequenza/2 + stepFreqRaffinato
nstepFrequenze = numpy.ceil((freqFin-freqIniz)/stepFreqRaffinato)+securbelt
#spindown
spindownIniz = headerSpindown[0]
stepSpindown = headerSpindown[1]
#nstepSpindown = 200
nstepSpindown = headerSpindown[2].astype(int)
print(stepSpindown)
# riarrangio gli array in modo che abbia i dati
# nel formato che voglio io
frequenze = tabella['frequenze'].values
frequenze = ((frequenze-freqIniz)/stepFreqRaffinato)-round(enhancement/2+0.001)
tempi = tabella['tempi'].values
tempi = tempi-epoca
tempi = ((tempi)*3600*24/stepFreqRaffinato)
#tempi = tempi - numpy.amin(tempi)+1
#tempi = tempi.astype(int)
print(minFreq, maxFreq, numpy.amin(tempi), numpy.amax(tempi))
print(maxFreq-minFreq, (numpy.amax(tempi)-numpy.amin(tempi))/1e11)
pesi = tabella['pesi'].values
#%reset_selective tabella
#nstepSpindown = 85
spindowns = numpy.arange(0, nstepSpindown)
spindowns = numpy.multiply(spindowns,stepSpindown)
spindowns = numpy.add(spindowns, spindownIniz)
# così ho i tre array delle tre grandezze, più i pesi e la fascia di sicurezza
#ora uso Tensorflow
securbeltTF = tf.constant(securbelt,dtype=tf.float32)
tempiTF = tf.constant(tempi,dtype=tf.float32)
pesiTF = tf.constant(pesi,dtype=tf.float32)
spindownsTF = tf.constant(spindowns, dtype=tf.float32)
frequenzeTF = tf.constant(frequenze, dtype=tf.float32)
nRows = tf.constant(nstepSpindown, dtype=tf.int32)
nColumns = nstepFrequenze
pesiTF = tf.reshape(pesiTF,(1,tf.size(pesi)))
pesiTF = pesiTF[0]
#Version con bincount, va leggermente più veloce su cpu ma
#attualmente non ha supporto gpu
def mapbincount(stepIesimo):
sdTimed = tf.multiply(spindowns[stepIesimo], tempi, name = "Tdotpert")
#sdTimed = tf.cast(sdTimed, dtype=tf.float32)
appoggio = tf.round(frequenze-sdTimed+securbelt/2, name = "appoggioperindici")
appoggio = tf.cast(appoggio, dtype=tf.int32)
valori = tf.bincount(appoggio,weights=pesi)
zeriDopo = tf.zeros([nColumns - tf.size(valori)], dtype=tf.float32)
riga = tf.concat([valori,zeriDopo],0, name = "rigadihough")
return riga
def mapnonVar(stepIesimo):
sdTimed = tf.multiply(spindownsTF[stepIesimo], tempiTF, name = "Tdotpert")
#sdTimed = tf.cast(sdTimed, dtype=tf.float32)
appoggio = tf.round(frequenzeTF-sdTimed+securbeltTF/2, name = "appoggioperindici")
appoggio = tf.cast(appoggio, dtype=tf.int32)
valori = tf.unsorted_segment_sum(pesiTF, appoggio, nColumns)
# zeriDopo = tf.zeros([nColumns - tf.size(valori)], dtype=tf.float32)
# riga = tf.concat([valori,zeriDopo],0, name = "rigadihough")
return valori
imagenonVar = tf.map_fn(mapnonVar, tf.range(0, nRows), dtype=tf.float32, parallel_iterations=8)
#sessione = tf.Session(config=tf.ConfigProto(log_device_placement=True))
sessione = tf.Session()
#run_options = tf.RunOptions(trace_level=tf.RunOptions.FULL_TRACE)
#run_metadata = tf.RunMetadata()
start = time.time()
#image = sessione.run(imagenonVar, options=run_options, run_metadata=run_metadata)
image = sessione.run(imagenonVar)
stop = time.time()
print(stop-start)
# Create the Timeline object, and write it to a json
#tl = timeline.Timeline(run_metadata.step_stats)
#ctf = tl.generate_chrome_trace_format()
#with open('timelinenonVar.json', 'w') as f:
# f.write(ctf)
nColumns = nColumns.astype(int)
semiLarghezza = numpy.round(enhancement/2+0.001).astype(int)
image[:,semiLarghezza*2:nColumns]=image[:,semiLarghezza*2:nColumns]-image[:,0:nColumns - semiLarghezza*2]
image = numpy.cumsum(image, axis = 1)
In [10]:
print(image.shape)
In [105]:
#FIN QUI SOLO HOUGH, ORA CALCOLO CANDS
#per ora numpy
numCand = 10
#if mode = 2
minDistance = enhancement*4
candidati = tf.zeros((9,20))
size = tf.shape(image, out_type=tf.float32)[1]
stepFreqRefTF = tf.constant(stepFreqRaffinato, dtype=tf.float32)
primaFreq = tf.constant(freqIniz-(securbelt/2)*stepFreqRaffinato, dtype=tf.float32)
freqniu = tf.range(size, delta = 1.0, dtype=tf.float32)
freqniu = freqniu*stepFreqRefTF+primaFreq
freqIniziale = tf.constant(scipy.io.loadmat(percorsoFile)['basic_info'][0,0]['frin'][0,0],dtype=tf.float32)
freqFinale = tf.constant(scipy.io.loadmat(percorsoFile)['basic_info'][0,0]['frfi'][0,0], dtype=tf.float32)
#QUI ANALOGO FUNZIONE CUT GD2
#%time indexInizialewh = numpy.where(freqniu>freqIniziale)[0][0]
#%time indexFinalewh = numpy.where(freqniu>freqFinale)[0][0]
indexIniziale = (freqIniziale-primaFreq)/stepFreqRefTF
indexIniziale = tf.cast(indexIniziale, dtype=tf.int32)
indexFinale = (freqFinale-primaFreq)/stepFreqRaffinato+1
indexFinale = tf.cast(indexFinale, dtype=tf.int32)
imageTF = tf.constant(image)
imageCand = imageTF[:,indexIniziale:indexFinale]
maxPerColumn = tf.reduce_max(imageCand, axis = 0)
rigaMax = tf.argmax(imageCand, axis = 0)
########################################################
stepFreqNew = tf.size(maxPerColumn)/numCand
stepFreqNew = tf.cast(stepFreqNew, dtype=tf.float32)
nMax = tf.size(maxPerColumn)
nMax = tf.cast(nMax, tf.float32)
indiciFreq = tf.range(1.0, nMax, delta=stepFreqNew)#, dtype=tf.int32)
indiciFreq = tf.concat([indiciFreq, tf.reshape(nMax+1, [1,1])[0]], 0)
indiciFreq = tf.cast(tf.round(indiciFreq), dtype=tf.int32)
print(sessione.run(indiciFreq).size)
def stats(tensor):
def get_median(v):
v = tf.reshape(v, [-1])
m = tf.shape(v)[0]//2
med = tf.nn.top_k(v, m).values[m-1]
return med
mediana = get_median(tensor)
sigmana= get_median(tf.abs(tensor-mediana))/0.6745 #cioè in unità dì sigma
#nb piccole differenze per float 32 invece di 64
stat = tf.stack([mediana,sigmana])
return stat
statTot = stats(imageCand)
#medianaTot = stats(imageCand)[0]
àsigmanaTot = stats(imageCand)[1]
medianaPerCand = tf.zeros(numCand)
sigmanaPerCand = tf.zeros(numCand)
iniziali = tf.concat([[indiciFreq[0]],indiciFreq[0:numCand-2],[indiciFreq[tf.size(indiciFreq)-4]]],0)
finali = tf.concat([[indiciFreq[2]-1],indiciFreq[3:numCand]-1,[indiciFreq[tf.size(indiciFreq)-2]-1]],0)
print(sessione.run(iniziali))
print(sessione.run(finali))
def statsPerCand(i):
stat = stats(maxPerColumn[iniziali[i]:finali[i]])#[0]
return stat
statPerCand = tf.map_fn(fn=statsPerCand,elems=tf.range(numCand-1),dtype=tf.float32)
sessione.run(statPerCand)
Out[105]:
In [ ]:
In [ ]:
for i in numpy.arange(1,numCand-1):
inizio = indiciFreq[i-1]
fine = indiciFreq[i+2]-1
medianaPerCand[i] = stats(maxPerColumn[inizio:fine],0.01)[0]
sigmanaPerCand[i] = stats(maxPerColumn[inizio:fine],0.01)[1]
i = 1
inizio = indiciFreq[i-1]
fine = indiciFreq[i+1]-1
medianaPerCand[0] = stats(maxPerColumn[inizio:fine],0.01)[0]
sigmanaPerCand[0] = stats(maxPerColumn[inizio:fine],0.01)[1]
i = indiciFreq.size-2
inizio = indiciFreq[i-2]
fine = indiciFreq[i]-1
medianaPerCand[numCand-1] = stats(maxPerColumn[inizio:fine],0.01)[0]
sigmanaPerCand[numCand-1] = stats(maxPerColumn[inizio:fine],0.01)[1]