In [2]:
import tensorflow as tf
import numpy
import scipy.io
import time
sessione = tf.Session()
# CARICO DATI
#tFft = 8192
tFft = 4096
tObs = 1/5 #mesi
tObs = tObs*30*24*60*60
nPunti = 2
cands = 100
primaFreq = 1/tFft
#percorsoDati = "dati/dati9mesi108HWI.mat"
#percorsoQuad = "quadHWI108.mat"
#percorsoPatch = "quadHWI108Ecl.mat"
#percorsoDati = "dati/dati9mesi108HWI.mat"
#percorsoQuad = ("quad%dLIL.mat" % tFft)
#percorsoPatch = ("quad%dEclNew.mat" % tFft)
percorsoDati = "/home/protoss/Documenti/TESI/wn100bkp/data/amplitude_4.mat"
#index = 9002 #signal
index = 4 #signal meno bello
#index = 5001
#index = 0
struttura1 = scipy.io.loadmat(percorsoDati)['H']
peakmap1 = struttura1[index].copy()
del struttura1
peakmapTOT = peakmap1
peaklinee = numpy.zeros(peakmapTOT.shape)
In [3]:
#spindowns
spindownMin = -10e-10
spindownMax = 9e-10
nstepSpindown = 100
stepSpindown = (spindownMax-spindownMin)/nstepSpindown
#stepSpindown = stepFreq/tObs
#nstepSpindown = numpy.round((spindownMax-spindownMin)/stepSpindown).astype(numpy.int32)
spindowns = numpy.arange(0, nstepSpindown)
spindowns = numpy.multiply(spindowns,stepSpindown)
spindowns = numpy.add(spindowns, spindownMin)
peaklinee = numpy.zeros(peakmapTOT.shape)
peaklinee[200] = 1
colonne = numpy.arange(peaklinee.shape[1])
print(colonne.size)
righe = numpy.round(numpy.arange(105,117,(117-105)/128)).astype(int)
print(righe.size)
peaklinee[righe,colonne] = 1
peaklinee = numpy.flip(peaklinee,0)
#peaklinee[:,10] = 1
In [5]:
from matplotlib import pyplot
import matplotlib as mpl
%matplotlib notebook
mpl.rcParams['font.size'] = 21
pyplot.figure(figsize=(10, 8))
#a = pyplot.scatter(numpy.arange(tempiUnici.size),tempiUnici/30, s = 0.5)
a = pyplot.imshow(peaklinee, cmap = 'gray_r',origin = 'lower', interpolation = 'none',aspect = 0.9)
pyplot.xlabel('$x$',fontsize = 21)
pyplot.ylabel('$y$', fontsize = 21)
pyplot.title('Input space', fontsize = 24)
pyplot.savefig('simpeak.pdf', format='pdf')
In [6]:
import tensorflow as tf
import numpy
import scipy.io
import time
sessione = tf.Session()
# CARICO DATI
#tFft = 8192
percorsoDati = "/home/protoss/Documenti/TESI/wn100bkp/data/amplitude_4.mat"
#index = 9002 #signal
index = 4 #signal meno bello
#index = 5001
#index = 0
struttura1 = scipy.io.loadmat(percorsoDati)['H']
#peakmap1 = struttura1[index].copy()
#del struttura1
struttura2 = scipy.io.loadmat(percorsoDati)['L']
#peakmap2 = struttura2[index].copy()
#del struttura2
struttura3 = scipy.io.loadmat(percorsoDati)['V']
#peakmap3 = struttura3[index].copy()
#del struttura3
strutturaTOT = struttura1 + struttura2 + struttura3
###############################
def noncorr():
freqMin = numpy.amin(frequenze)
freqIniz = freqMin - stepFreq/2 -stepFreqRaffinato
freqNonCor = (frequenze -freqIniz)/stepFreqRaffinato-round(enhancement/2+0.001)
#freqNonCor = tf.constant(freqNonCor, dtype = tf.float32)
return freqNonCor, freqIniz
def inDaHough(freqHM):
#calcola la hough per ogni step di spindown
def houghizza(stepIesimo):
sdTimed = tf.multiply(spindownsTF[stepIesimo], tempiHM, name = "Tdotpert")
#sdTimed = tf.cast(sdTimed, dtype=tf.float32)
appoggio = tf.round(freqHM-sdTimed+securbelt/2, name = "appoggioperindici")
appoggio = tf.cast(appoggio, dtype=tf.int32)
valorisx = tf.unsorted_segment_sum(pesiHM, appoggio, nColumns)
return valorisx
houghDiff = tf.map_fn(houghizza, tf.range(0, nRows), dtype=tf.float32, parallel_iterations=8)
def sliceInt():
#faccio integrazione finale (vecchia versione senza conv)
semiLarghezza = tf.round(enhancement/2+0.001)
semiLarghezza = tf.cast(semiLarghezza, tf.int32)
houghInt = houghDiff[:,enhancement:nColumns]-houghDiff[:,0:nColumns - enhancement]
houghInt = tf.concat([houghDiff[:,0:enhancement],houghInt],1)
return houghInt
hough = sliceInt()
#hough = convInt()
houghinal = tf.cumsum(hough, axis = 1)
return houghinal
#return houghDiff
tFft = 4096
tObs = 1/5 #mesi
tObs = tObs*30*24*60*60
#nPunti = 2
cands = 32
primaFreq = 1/tFft
#maximumCand = numpy.array()
peakmapTOT = strutturaTOT[index]
nonzeri = numpy.nonzero(peakmapTOT)
peakmapTOT[nonzeri] = 1
In [7]:
#sparsa = numpy.nonzero(peakmapTOT)
sparsa = numpy.nonzero(peaklinee)
frequenze,tempi = sparsa
tempi = tempi+1
frequenze = frequenze / tFft + 1
pesi = numpy.ones(sparsa[0].size)
#nb: picchi ha 0-tempi
# 1-frequenze
# 2-pesi
#headers vari
securbelt = 2000
#securbelt = 4000*3
#frequenze
stepFreq = 1/tFft
enhancement = 10
#enhancement = 1
stepFreqRaffinato = stepFreq/enhancement
nstepsFreq = securbelt+(numpy.amax(frequenze)-numpy.amin(frequenze) + stepFreq + 2*stepFreqRaffinato)/stepFreqRaffinato
epoca = 55
#spindowns
spindownMin = -9e-10
spindownMax = 9e-10
nstepSpindown = 100
stepSpindown = (spindownMax-spindownMin)/nstepSpindown
spindowns = numpy.arange(0, nstepSpindown)
spindowns = numpy.multiply(spindowns,stepSpindown)
spindowns = numpy.add(spindowns, spindownMin)
print(numpy.where(spindowns == numpy.amin(numpy.absolute(spindowns))))
#da qui si usa tensorflow
#definisco tutte le costanti necessarie
tempiTF = tf.constant(tempi, dtype=tf.float64)
pesiTF = tf.constant(pesi,dtype=tf.float32)
spindownsTF = tf.constant(spindowns, dtype=tf.float32)
tempiHM = tempiTF-epoca
tempiHM = ((tempiHM)*3600*24/stepFreqRaffinato)
tempiHM = tf.cast(tempiHM, tf.float32)
pesiHM = tf.reshape(pesiTF,(1,tf.size(pesiTF)))
pesiHM = pesiHM[0]
nRows = tf.constant(nstepSpindown, dtype=tf.int32)
nColumns = tf.cast(nstepsFreq, dtype=tf.int32)
freq, freqIn = noncorr()
freqTF = tf.constant(freq, dtype = tf.float32)
houghmap = inDaHough(freqTF)
hough = sessione.run(houghmap)
In [10]:
from matplotlib import pyplot
import matplotlib as mpl
%matplotlib notebook
mpl.rcParams['font.size'] = 22
pyplot.figure(figsize=(12, 8))
a = pyplot.imshow(hough[:,400:-50], cmap = 'gray_r', origin = 'lower', interpolation = 'none', aspect = 30)
labelxtick = [0,50,100,150,200,250]
posxTick = [0,500,1000,1500,2000,2500]
pyplot.xticks(posxTick,labelxtick)
sdveri = (-0.1953/15)*(numpy.arange(101)-50)
print(sdveri)
#labelytick = numpy.array([-50, -25, 0, 25, 50])
sdlabel = numpy.array([0,17,33,50,66,83,100])
labelytick = sdveri[sdlabel]
labelytick = [-0.65, -0.43, -0.21, 0, 0.21, 0.43, 0.65]
posytick = sdlabel#numpy.arange(0,101,25)#[ 0, 25 ,50, 75, 100]
pyplot.yticks(posytick, labelytick)
pyplot.colorbar(shrink = 0.5,aspect = 15)
pyplot.xlabel('$q$',fontsize = 22)
pyplot.ylabel('$m$',fontsize = 22)
pyplot.title('Parameters space', fontsize = 25)
pyplot.tight_layout()
pyplot.savefig('simphough.pdf', format='pdf')
In [8]:
colonne = numpy.arange(peakmapTOT.shape[1])
print(colonne.size)
righe = numpy.round(numpy.arange(135,151,(151-135)/128)).astype(int)
righe = numpy.flip(righe,0)
print(righe.size)
peakmapTOT[righe,colonne] = 1
In [9]:
from matplotlib import pyplot
%matplotlib notebook
pyplot.figure(figsize=(12, 8))
#peaklinee = numpy.flip(peaklinee,0)
nonzeri = numpy.nonzero(peakmapTOT)
peakmapTOT[nonzeri] = 1
a = pyplot.imshow(peakmapTOT, cmap = 'binary', origin = 'lower', interpolation = 'none', aspect = 0.9)
#sdveri = (-0.1953/15)*(numpy.arange(101)-50)
#print(sdveri)
#labelytick = numpy.array([-50, -25, 0, 25, 50])
#sdlabel = numpy.array([0,17,33,50,66,83,100])
#labelytick = sdveri[sdlabel]
#labelytick = [-0.65, -0.43, -0.21, 0, 0.21, 0.43, 0.65]
labelxtick = numpy.array([0,20.0,40.0,60.0,80.0,100.0,120.0])*8192/3600
print(labelxtick)
labelytick = numpy.array([0,50.0,100.0,150.0,200.0,250.0])
labelytick = labelytick/8192+80
print(labelytick)
pyplot.xlabel('$t$ $(hours)$')
labelxtick = [0,23,46,69,92,115,138]
posxTick = [0,20,40,60,80,100,120]
pyplot.xticks(posxTick,labelxtick)
pyplot.ylabel('$\\nu$ ($\\times 10^{-3} + 80 Hz$)' )
labelytick = [0,6,12,18,24,30]
posyTick = [0,50,100,150,200,250]
pyplot.yticks(posyTick,labelytick)
pyplot.title('Peakmap')
pyplot.show()
In [11]:
import tensorflow as tf
import numpy
import scipy.io
import time
sessione = tf.Session()
# CARICO DATI
#tFft = 8192
percorsoDati = "/home/protoss/Documenti/TESI/wn100bkp/data/amplitude_4.mat"
#index = 9002 #signal
index = 4 #signal meno bello
#index = 5001
#index = 0
struttura1 = scipy.io.loadmat(percorsoDati)['H']
#peakmap1 = struttura1[index].copy()
#del struttura1
struttura2 = scipy.io.loadmat(percorsoDati)['L']
#peakmap2 = struttura2[index].copy()
#del struttura2
struttura3 = scipy.io.loadmat(percorsoDati)['V']
#peakmap3 = struttura3[index].copy()
#del struttura3
strutturaTOT = struttura1 + struttura2 + struttura3
###############################
def noncorr():
freqMin = numpy.amin(frequenze)
freqIniz = freqMin - stepFreq/2 -stepFreqRaffinato
freqNonCor = (frequenze -freqIniz)/stepFreqRaffinato-round(enhancement/2+0.001)
#freqNonCor = tf.constant(freqNonCor, dtype = tf.float32)
return freqNonCor, freqIniz
def inDaHough(freqHM):
#calcola la hough per ogni step di spindown
def houghizza(stepIesimo):
sdTimed = tf.multiply(spindownsTF[stepIesimo], tempiHM, name = "Tdotpert")
#sdTimed = tf.cast(sdTimed, dtype=tf.float32)
appoggio = tf.round(freqHM-sdTimed+securbelt/2, name = "appoggioperindici")
appoggio = tf.cast(appoggio, dtype=tf.int32)
valorisx = tf.unsorted_segment_sum(pesiHM, appoggio, nColumns)
return valorisx
houghDiff = tf.map_fn(houghizza, tf.range(0, nRows), dtype=tf.float32, parallel_iterations=8)
def sliceInt():
#faccio integrazione finale (vecchia versione senza conv)
semiLarghezza = tf.round(enhancement/2+0.001)
semiLarghezza = tf.cast(semiLarghezza, tf.int32)
houghInt = houghDiff[:,enhancement:nColumns]-houghDiff[:,0:nColumns - enhancement]
houghInt = tf.concat([houghDiff[:,0:enhancement],houghInt],1)
return houghInt
hough = sliceInt()
#hough = convInt()
houghinal = tf.cumsum(hough, axis = 1)
return houghinal
#return houghDiff
tFft = 4096
tObs = 1/5 #mesi
tObs = tObs*30*24*60*60
#nPunti = 2
cands = 32
primaFreq = 1/tFft
#maximumCand = numpy.array()
nIndexes = 10
#nIndexes = strutturaTOT.shape[0]
allMaxCands = numpy.zeros((nIndexes,4,5))
peakmapTOT = strutturaTOT[index]
#del struttura1
#del struttura2
#del struttura3
colonne = numpy.arange(peakmapTOT.shape[1])
print(colonne.size)
righe = numpy.round(numpy.arange(135,151,(151-135)/128)).astype(int)
righe = numpy.flip(righe,0)
print(righe.size)
peakmapTOT[righe,colonne] = 1
sparsa = numpy.nonzero(peakmapTOT)
frequenze,tempi = sparsa
tempi = tempi+1
frequenze = frequenze / tFft + 1
pesi = numpy.ones(sparsa[0].size)
#nb: picchi ha 0-tempi
# 1-frequenze
# 2-pesi
#headers vari
securbelt = 2000
#securbelt = 4000*3
#frequenze
stepFreq = 1/tFft
enhancement = 10
#enhancement = 1
stepFreqRaffinato = stepFreq/enhancement
nstepsFreq = securbelt+(numpy.amax(frequenze)-numpy.amin(frequenze) + stepFreq + 2*stepFreqRaffinato)/stepFreqRaffinato
epoca = 55
#spindowns
spindownMin = -10e-10
spindownMax = 9e-10
nstepSpindown = 90
stepSpindown = (spindownMax-spindownMin)/nstepSpindown
spindowns = numpy.arange(0, nstepSpindown)
spindowns = numpy.multiply(spindowns,stepSpindown)
spindowns = numpy.add(spindowns, spindownMin)
#da qui si usa tensorflow
#definisco tutte le costanti necessarie
tempiTF = tf.constant(tempi, dtype=tf.float64)
pesiTF = tf.constant(pesi,dtype=tf.float32)
spindownsTF = tf.constant(spindowns, dtype=tf.float32)
tempiHM = tempiTF-epoca
tempiHM = ((tempiHM)*3600*24/stepFreqRaffinato)
tempiHM = tf.cast(tempiHM, tf.float32)
pesiHM = tf.reshape(pesiTF,(1,tf.size(pesiTF)))
pesiHM = pesiHM[0]
nRows = tf.constant(nstepSpindown, dtype=tf.int32)
nColumns = tf.cast(nstepsFreq, dtype=tf.int32)
freq, freqIn = noncorr()
freqTF = tf.constant(freq, dtype = tf.float32)
houghmap = inDaHough(freqTF)
hough = sessione.run(houghmap)
from matplotlib import pyplot
%matplotlib notebook
pyplot.figure(figsize=(12, 8))
a = pyplot.imshow(hough[:,1000:-1000], cmap = 'gray_r', origin = 'lower', interpolation = 'none', aspect = 30)
labelxtick = [0,6,12,18,24,30]
posxTick = [0,500,1000,1500,2000,2500]
#pyplot.xticks(posxTick,labelxtick)
labelytick = numpy.arange(-8,8.1,2)
print(labelytick)
posyTick = [7,17,27,37,47, 57,67,77,87]
#pyplot.yticks(posytick,labelytick)
pyplot.xlabel('$\\nu_0$ ($\\times 10^{-3} + 80 Hz$)')
pyplot.ylabel('$\dot{\\nu} \;(10^{-10}\:Hz/s)$')
pyplot.xticks(posxTick,labelxtick)
pyplot.yticks(posyTick,labelytick)
pyplot.colorbar(shrink = 0.5,aspect = 15)
pyplot.title('Frequency-Hough map')
pyplot.show()
In [35]:
zeri = numpy.where(hough[:,1000:-1000] == 0)
zeri
Out[35]:
In [59]:
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import axes3d, Axes3D #<-- Note the capitalization!
%matplotlib notebook
indici = numpy.nonzero(hough[:,1000:-1000])
righe = indici[0]
colonne = indici[1]
righe = numpy.concatenate([righe,[86,86]])
colonne = numpy.concatenate([colonne,[2529,2530]])
valori = numpy.ravel(hough[:,1000:-1000])
print(righe.size, colonne.size, valori.size)
fig = plt.figure()
ax = Axes3D(fig) #<-- Note the difference from your original code...
# Plot the values
ax.scatter(righe, colonne, valori, c = valori, marker='o', alpha = 1)
Out[59]:
In [26]:
from scipy import sparse
sparsh = sparse.csr_matrix(hough[:,1000:-1000])
In [30]:
print(sparsh.shape)
In [ ]: