際際滷

際際滷Share a Scribd company logo

Test 03a 08.05.2015

Corrado Pecora
Corrado Pecora

Continuazione di Test 3 del 21.04.2015

1 of 24
Downloaded 29 times
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 1 di 24 TEST 3a Acquisitore dati - Progetto DOLFRANG - Prototipo 5 open source Acquisizione del 20/04/2015 Continuazione di Test 3 del 21.04.2015 Si prosegue lillustrazione dei risultati ottenuti elaborando i dati dellacquisizione. Il testo precedente lasciava in sospeso la questione dellequalizzazione del segnale. E da qui si riprende, riportando lultima parte del testo precedente. Per quanto riguarda i risultati dellelaborazione relativa allo spettro di risposta sui singoli canali, occorre lequalizzazione del segnale per tenere conto della risposta del sensore in relazione alla frequenza. Lo spettro di ciascun canale senza equalizzazione 竪 riportato nel seguito, relativamente al caso con gain 2.100 (24k ohm). In tali risultati la parte al di sotto di 7Hz richiede equalizzazione. Note di elaborazione: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti-triggering on raw signal, n属 finestre 134
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 2 di 24 Ricordando che la curva di risposta dei sensori geofoni 竪 si nota il decadimento del segnale per valori inferiori a 7Hz. Quindi, per analizzare i risultati con riferimento al singolo canale, occorre procedere ad una correzione dei valori del segnale acquisito per tenere conto della dinamica del geofono utilizzato. Come si vedr nel seguito, questa operazione non 竪 necessaria per le valutazioni H/V, in quanto si tratta di rapporto tra quantit soggette al medesimo fattore di scala in termini di frequenza, quindi il risultato finale 竪 sovrapponibile. Per equalizzare il segnale sono state svolte le seguenti operazioni: 1) interpolazione numerica della curva di risposta del geofono, in modo da avere un modello matematico applicabile nei passi successivi; 2) calcolo della Trasformata di Fourier del segnale acquisito;
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 3 di 24 3) applicazione del modello matematico numerico della curva di risposta del geofono, di cui al punto (1), nel dominio delle frequenze del segnale in elaborazione; 4) calcolo della Antitrasformata di Fourier per ricostruire il segnale corretto nel dominio del tempo. Per loperazione (1) 竪 stato utilizzato un foglio elettronico, facilmente realizzabile o con software free quali OpenOffice, LibreOffice, ecc. o commerciale quale Excel. Per le operazioni descritte ai punti (2), (3) e (4) 竪 stato utilizzato un ambiente matematico denominato Scilab. Dal sito http://www.scilab.org/ : Scilab 竪 un pacchetto di programmi gratuiti per la computazione numerica sviluppati dallo INRIA e dallo ENPC inFrancia, poi da Scilab Consortium in seno alla Fondazione Digiteo. Oggi Scilab 竪 sviluppato da Scilab Enterprises. Scilab attualmente utilizza la licenza CeCILL v2, che 竪 compatibile con la licenza GNU GPL v2 della Free Software Foundation. La versione utilizzata per le analisi 竪 la 5.5.2 sia in ambiente a 32 bit sia in quello a 64 bit (s.o. Windows 7 pro). Tale scelta 竪 stata fatta in applicazione della Mission del sistema Theremino, ossia: Ricongiungere il mondo digitale col mondo reale e concreto. Con tecnologie semplici e facilmente riproducibili. E con il massimo rispetto per lambiente e per tutti gli esseri viventi. In questo caso, mediante un acquisitore autocostruito, 竪 stato campionato un segnale relativo a un microtremore. La sua elaborazione 竪 stata finora condotta in un ambiente software specializzato, quale Geopsy. Adesso si mostra come iniziare ad operare in un ambiente software pi湛 generico, quale Scilab, rendendo ancora pi湛 evidente la possibilit di sfruttare a pieno la potenzialit dei computer, volendo anche con alternative free a Scilab quali Octave, SciPy, python(x,y), ecc. o commerciali quali Matlab, Mathematica, ecc.. In appendice viene riportato il listato delle istruzioni utilizzate in ambiente Scilab. Nel seguito si illustrano i risultati ottenuti. Nella figura seguente si riporta (1) linterpolazione della curva di risposta del geofono
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 4 di 24 DigitalizzazionedellacurvadirispostadelgeofonoSensheSS-4.5N frequenzarispostav/(m/s)decadimentodecadimento% logXlogYXY(28,8-Y)(28,8-Y)/28,8% 0,0000,1501,0001,4127,4395,10 0,5151,1773,27315,0313,8147,88 0,5791,2803,79319,059,7933,93 0,6171,3304,14021,387,4625,87 0,6561,3704,52923,445,4018,72 0,6921,4004,92025,123,7212,90 0,7431,4305,53426,921,926,67 0,7891,4606,15228,840,000,00 0,8461,4607,01528,840,000,00 1,7161,46052,00028,840,000,00 2,3051,460201,83728,840,000,00 2,7591,460574,11628,840,000,00 2,9651,460922,57128,840,000,00 3,0001,4601000,00028,840,000,00 Interpolazionematematicapernpunticonpolinomiodigradon-1 n属XY 11,0001,413Controllointerp.012345 23,27315,031yta0a1a2a3a4a5 33,79319,055a0221,890954619,0551,00E+003,79E+001,44E+015,46E+012,07E+021,22E+03 44,14021,380a1-192,357807821,3801,00E+004,14E+001,71E+017,10E+012,94E+021,22E+03 54,52923,442a264,8499648623,4421,00E+004,53E+002,05E+019,29E+014,21E+021,91E+03 64,92025,119a3-9,30523319125,1191,00E+004,92E+002,42E+011,19E+025,86E+022,88E+03 75,53426,915a40,49615233526,9151,00E+005,53E+003,06E+011,69E+029,38E+025,19E+03 86,15228,840a5-0,00092490828,8401,00E+006,15E+003,78E+012,33E+021,43E+038,81E+03 97,01528,840 1052,00028,840 11201,83728,840 12574,11628,840 13922,57128,840 141000,00028,840 Parametri Rettada(1)a(3) Rettada(8)a(14) Matricedivalutazionedelpolinomio 1 1 2 210 ++++=n nxa...xaxaay
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 5 di 24 ossia: (a) una retta nel piano log10-log10; (b) una interpolazione polinomiale del 5属 ordine per il tratto di raccordo; (c) un tratto di plateau con risposta costante del geofono (valore massimo 28.84 V/m/s). Per quanto riguarda i passi di calcolo (2), (3) e (4) nel seguito si illustrano i risultati dellelaborazione ottenuti in ambiente Scilab con riferimento al segnale acquisito con gain 2100 (24k Ohm). Passo di calcolo (2), la Trasformata di Fourier del segnale acquisito:
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 6 di 24 Passo di calcolo (3), lequalizzazione della Trasformata di Fourier con linterpolazione della curva di risposta del geofono
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 7 di 24
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 8 di 24 Passo di calcolo (4), la ricostruzione del segnale nel dominio del tempo
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 9 di 24 Infine un confronto tra un estratto del file .saf originario e di quello ottenuto a seguito dellequalizzazione, che consente di capire la compatibilit del nuovo file con Geopsy segnale originario segnale equalizzato
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 10 di 24 Utilizzando Geopsy si espongono i risultati dellelaborazione in termini di spettro di risposta dei singoli canali Segnale originario: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti- triggering on raw signal, n属 finestre 134 Segnale equalizzato: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti- triggering on raw signal, n属 finestre 122
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 11 di 24 Utilizzando la funzione Spectrum Rotate si ottiene Segnale originario: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti- triggering on raw signal, campionamento 200 punti, n属 finestre 41 Segnale equalizzato: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti- triggering on raw signal, campionamento 200 punti, n属 finestre 20 Nel caso specifico si osservi che lequalizzazione delle basse frequenze poco ha aggiunto alle analisi gi svolte in termini di individuazione delle frequenze proprie dominanti del manufatto.
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 12 di 24 Infine utilizzando la funzione H/V si ottiene Segnale originario: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti-triggering on raw signal, n属 finestre 24 Segnale equalizzato: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti-triggering on raw signal, n属 finestre 20 La sostanziale coincidenza dei risultati conferma quanto esposto in precedenza sulla non necessit di equalizzare un segnale per analisi H/V. Quindi adesso il processo di calcolo 竪 abbastanza completo sulla manipolazione del singolo segnale nel dominio delle frequenze. Nei successivi test si mostreranno altri aspetti di calcolo e comportamento applicati bagli edifici. Grazie per lattenzione.
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 13 di 24 APPENDICE 1 Evoluzione dellacquisitore USB Lacquisitore Hear Our Earth ha avuto una evoluzione: 1 per ogni canale 竪 possibile selezionare un diverso valore di gain, correlato al valore di resistenza posizionata tra il geofono e lamplificatore. Le resistenze selezionabili sono 0 Ohm, 150k Ohm e 24k Ohm. La selezione avviene mediante commutatore rotativo; 2 il computer per lacquisizione 竪 un tablet Asus M81C con sistema operativo Windows 8.1, che consente lacquisizione contemporanea di almeno 6 canali a 500Hz (non avendo altri acquisitori a disposizione per un testare un maggiore carico). Per agevolare limmissione dei dati si pu嘆 utilizzare una tastiera/mouse bluetooth Ipazzport- KP-810-19BTT. Ovviamente 竪 occorso un po di lavoro per le modifiche, ma questo fa parte della volont di sperimentare.
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 14 di 24 alcune fasi delle modifiche nota: il corto con 375 Ohm 竪 gi compreso nelle caratteristiche del geofono
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 15 di 24 APPENDICE 2 Equalizzazione del segnale Software: Scilab ver. 5.5.2 con Windows 7 pro a 64 bit Utilizzo: sulla console di Scilab vengono visualizzati i passi di calcolo, illustrati nel seguito. In caratteri di colore nero evidenziano quanto appare a video, con quelli di colore rosso e in corsivo si riporta un commento al singolo passo di calcolo. *** FASE 1 DI 8: SELEZIONE FILE DATI *** mediante una maschera di windows 竪 possibile selezionare il file dati con estensione .saf *** INIZIO *** *** FINE *** *** FASE 2 DI 8: LETTURA FILE DATI *** vengono letti i dati dal file *** INIZIO *** *** FINE *** *** FASE 3 DI 8: CALCOLO FFT *** calcolo della Fast Fourier Trasformcon la libreria fftw3, riferimenti in http://www.fftw.org/ *** INIZIO *** *** FINE *** *** FASE 4 DI 8: VISUALIZZAZIONE FFT *** grafici dei segnali relativi ai 3 canali Verticale, Nord e Est *** INIZIO *** *** FINE *** *** FASE 5 DI 8: EQUALIZZAZIONE *** applicazione del modello numerico della curva di decadimento del geofono *** INIZIO *** *** FINE *** *** FASE 6 DI 8: RICOSTRUZIONE DEL SEGNALE *** calcolo della Antitrasformata di Fourier al segnale equalizzato nel dominio delle frequenze *** INIZIO *** *** FINE *** *** FASE 7 DI 8: SCRITTURA SEGNALE CORRETTO *** scrittura su file (con nome modificato rispetto loriginale con laggiunta in coda di _EQ.saf) *** INIZIO *** *** FINE *** *** FASE 8 DI 8: CANCELLAZIONE VARIABILI *** rilascio della memoria utilizzata *** INIZIO *** *** FINE ***
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 16 di 24 Listato: // Equalizzazione di un segnale geofonico // Geofoni SS-4.5N della Senshe //---------------------------------- // Data: 08.05.2015 // Autore: Corrado Pecora mprintf('%s n',"*** INIZIO ELABORAZIONE ***") // Inizializzazione della memoria // Sistema Operativo: Windows 7 Pro 64 bit // Ram: 16 Gb stacksize(100000000) // 1 - SELEZIONE FILE DATI mprintf('%s n',"*** FASE 1 DI 8: SELEZIONE FILE DATI ***") // maschera di dialogo windows // selezione di un solo file SAF mprintf('%s n',"*** INIZIO ***") FileName=uigetfile(["*.saf"],"","Seleziona un file SAF",%f) mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 2 - LETTURA DATI DAL FILE mprintf('%s n',"*** FASE 2 DI 8: LETTURA FILE DATI ***") mprintf('%s n',"*** INIZIO ***") // apertura file testo in lettura fd_r=mopen(FileName,'rt') // lettura delle prime 9 righe for i=1:9 riga=mgetl(fd_r,1) end // lettura della frequenza di campionamento // (sulla stessa riga: stringa %s, stringa %s, intero %g) riga=mfscanf(fd_r,'%s %s %g') SAMP_FREQ=riga(:,3) // lettura del numero di dati // (sulla stessa riga: stringa %s, stringa %s, intero %g) riga=mfscanf(fd_r,'%s %s %g') NDAT=riga(:,3)
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 17 di 24 // lettura di 9 righe (una per EndOfLine della precedente) for i=1:9 riga=mgetl(fd_r,1) end // lettura di terne di dati in NDAT righe relativi a velocit geofoniche segnale(1:NDAT,1:3)=0.0 for i=1:NDAT segnale(i,1:3)=mfscanf(fd_r,'%g %g %g') end // chiusura file testo mclose(fd_r) // asse dei tempi del segnale acquisito // ascisse equispaziate di dt = 1/SAMP_FREQ dt=1.0/SAMP_FREQ tempo = 0.0:dt:(NDAT-1)*dt mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 3 - FFT mprintf('%s n',"*** FASE 3 DI 8: CALCOLO FFT ***") mprintf('%s n',"*** INIZIO ***") // suddivisione del segnale segnale_V=segnale(1:NDAT,1) segnale_N=segnale(1:NDAT,2) segnale_E=segnale(1:NDAT,3) // calcolo fft_segnale_V=fft(segnale_V,-1,"symmetric") fft_segnale_N=fft(segnale_N,-1,"symmetric") fft_segnale_E=fft(segnale_E,-1,"symmetric") mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 4 - VISUALIZZAZIONE DELL'ELABORAZIONE mprintf('%s n',"*** FASE 4 DI 8: VISUALIZZAZIONE FFT ***") mprintf('%s n',"*** INIZIO ***") // essendo "segnale" un vettore di numeri reali // la fft sar coniugata simmetrica. Quindi verranno visualizzati solo
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 18 di 24 // i primi NDAT/2 punti. // frequenze associate // VEDI http://usingscilab.blogspot.it/2009/08/fft-specifying-frequencies.html // http://www.mathworks.com/matlabcentral/answers/36430-plotting-the-frequency-spectrum dF=SAMP_FREQ/NDAT frequenza=(0:NDAT-1)*dF NrDati=size(frequenza,'*') freq_plot=frequenza(1:NrDati/2) // i 3 assi distinti VVV=fft_segnale_V(1:NrDati/2); //asse Verticale NNN=fft_segnale_N(1:NrDati/2); //asse Nord EEE=fft_segnale_E(1:NrDati/2); //asse Est // grafici (nota: frequenza' 竪 la trasposta di frequenza) scf(1) clf(1) subplot(3,1,1) plot(freq_plot',real(VVV)) xtitle("Asse Verticale - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(VVV)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(VVV)) xtitle("Valore assoluto") scf(2) clf(2) subplot(3,1,1) plot(freq_plot',real(NNN)) xtitle("Asse Nord - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(NNN)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(NNN)) xtitle("Valore assoluto")
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 19 di 24 scf(3) clf(3) subplot(3,1,1) plot(freq_plot',real(EEE)) xtitle("Asse Est - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(EEE)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(EEE)) xtitle("Valore assoluto") mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 5 - EQUALIZZAZIONE mprintf('%s n',"*** FASE 5 DI 8: EQUALIZZAZIONE ***") mprintf('%s n',"*** INIZIO ***") // Nota: la curva di equalizzazione 竪 stata ricavata mediante // digitalizzazione del grafico della risposta in frequenza // parametri della interpolazione // 5.1) limiti in frequenza dei 3 campi: rettilineo, polinomiale, // rettilineo orizzontale // N.B.: la retta 竪 nel piano log10-log10 come la rappresentazione della // curva di risposta in frequenza del geofono freq_0=1.00 freq_1=3.79315 freq_2=6.09 // minima frequenza di equalizzazione freq_min=0.3 //max valore di risposta in frequenza ymax=28.84 // 5.2) punti della retta del campo (1) x0=log10(freq_0) x1=log10(freq_1) y0=log10(1.413) y1=log10(19.055)
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 20 di 24 m=(y1-y0)/(x1-x0) // 5.3) parametri della interpolazione polinomiale del 5属 ordine a0=221.8909546 a1=-192.3578078 a2=64.84996486 a3=-9.305233191 a4=0.496152335 a5=-0.000924908 // 5.4) correzione della FFT per ogni singolo canale fft1_segnale_V(1:NrDati)=fft_segnale_V(1:NrDati) fft1_segnale_N(1:NrDati)=fft_segnale_N(1:NrDati) fft1_segnale_E(1:NrDati)=fft_segnale_E(1:NrDati) // ciclo di correzione dei valori i=2; //per i=1 frequenza(i)=0.0 freq_=frequenza(i) while freq_< freq_2 if freq_<= freq_min then // per questioni numeriche si equalizza la fft a partire da una certa // frequenza Coeff=1.0 else if freq_<=freq_1 then // tratto rettilineo Coeff=ymax/10.^((log10(freq_)-x0)*m+y0) else // interpolazione polinomiale Coeff=ymax/(a0+a1*freq_+a2*freq_^2+a3*freq_^3+a4*freq_^4+a5*freq_^5) end end fft1_segnale_V(i)=Coeff*fft_segnale_V(i) fft1_segnale_N(i)=Coeff*fft_segnale_N(i) fft1_segnale_E(i)=Coeff*fft_segnale_E(i) i=i+1 freq_=frequenza(i) end // i 3 assi distinti
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 21 di 24 VVV1=fft1_segnale_V(1:NrDati/2); //asse Verticale NNN1=fft1_segnale_N(1:NrDati/2); //asse Nord EEE1=fft1_segnale_E(1:NrDati/2); //asse Est // 5.5) grafici scf(4) clf(4) subplot(3,1,1) plot(freq_plot',real(VVV1)) xtitle("Asse Verticale - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(VVV1)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(VVV1)) xtitle("Valore assoluto") scf(5) clf(5) subplot(3,1,1) plot(freq_plot',real(NNN1)) xtitle("Asse Nord - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(NNN1)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(NNN1)) xtitle("Valore assoluto") scf(6) clf(6) subplot(3,1,1) plot(freq_plot',real(EEE1)) xtitle("Asse Est - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(EEE1)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(EEE1)) xtitle("Valore assoluto")
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 22 di 24 mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 6 - RICOSTRUZIONE DEL SEGNALE mprintf('%s n',"*** FASE 6 DI 8: RICOSTRUZIONE DEL SEGNALE ***") mprintf('%s n',"*** INIZIO ***") // Calcolo della trasformata inversa di Fourier sulla base della // trasformata corretta segnale1_V=fft(fft1_segnale_V,1) //asse Verticale segnale1_N=fft(fft1_segnale_N,1) //asse Nord segnale1_E=fft(fft1_segnale_E,1) //asse Est // Grafici scf(7) clf(7) subplot(2,1,1) plot(tempo',segnale_V) xtitle("Asse Verticale - Segnale originario") subplot(2,1,2) plot(tempo',segnale1_V) xtitle("Asse Verticale - Segnale corretto") scf(8) clf(8) subplot(2,1,1) plot(tempo',segnale_N) xtitle("Asse Nord - Segnale originario") subplot(2,1,2) plot(tempo',segnale1_N) xtitle("Asse Nord - Segnale corretto") scf(9) clf(9) subplot(2,1,1) plot(tempo',segnale_E) xtitle("Asse Est - Segnale originario") subplot(2,1,2) plot(tempo',segnale1_E) xtitle("Asse Est - Segnale corretto")
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 23 di 24 mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 7 - SCRITTURA DEL SEGNALE CORRETTO mprintf('%s n',"*** FASE 7 DI 8: SCRITTURA SEGNALE CORRETTO ***") mprintf('%s n',"*** INIZIO ***") // inserimento di"_EQ" alla fine del nome del file di input FileNameOut=strsubst(FileName,".saf","_EQ.saf") // apertura file input in lettura fd_r=mopen(FileName,'rt') // apertura del file output in scrittura fd_w=mopen(FileNameOut,'wt') // lettura e scrittura delle prime 19 righe for i=1:19 riga=mgetl(fd_r,1) mputl(riga,fd_w) end // chiusura del file di input mclose(fd_r) // scrittura di terne di valori per ogni riga for i=1:NDAT riga=sprintf('%8.3f %8.3f %8.3f ',segnale1_V(i),segnale1_N(i),segnale1_E(i)) mputl(riga,fd_w) end // chiusura del file di output mclose(fd_w) mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 8 - CANCELLAZIONE DELLE VARIABILI // per rilascio della memoria mprintf('%s n',"*** FASE 8 DI 8: CANCELLAZIONE VARIABILI ***") mprintf('%s n',"*** INIZIO ***") clear
Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 24 di 24 mprintf('%s n',"*** FINE ***") mprintf('%s n',"") mprintf('%s n',"*** FINE ELABORAZIONE ***") mprintf('%s n',"")

More Related Content

Test 03a 08.05.2015

  • 1. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 1 di 24 TEST 3a Acquisitore dati - Progetto DOLFRANG - Prototipo 5 open source Acquisizione del 20/04/2015 Continuazione di Test 3 del 21.04.2015 Si prosegue lillustrazione dei risultati ottenuti elaborando i dati dellacquisizione. Il testo precedente lasciava in sospeso la questione dellequalizzazione del segnale. E da qui si riprende, riportando lultima parte del testo precedente. Per quanto riguarda i risultati dellelaborazione relativa allo spettro di risposta sui singoli canali, occorre lequalizzazione del segnale per tenere conto della risposta del sensore in relazione alla frequenza. Lo spettro di ciascun canale senza equalizzazione 竪 riportato nel seguito, relativamente al caso con gain 2.100 (24k ohm). In tali risultati la parte al di sotto di 7Hz richiede equalizzazione. Note di elaborazione: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti-triggering on raw signal, n属 finestre 134
  • 2. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 2 di 24 Ricordando che la curva di risposta dei sensori geofoni 竪 si nota il decadimento del segnale per valori inferiori a 7Hz. Quindi, per analizzare i risultati con riferimento al singolo canale, occorre procedere ad una correzione dei valori del segnale acquisito per tenere conto della dinamica del geofono utilizzato. Come si vedr nel seguito, questa operazione non 竪 necessaria per le valutazioni H/V, in quanto si tratta di rapporto tra quantit soggette al medesimo fattore di scala in termini di frequenza, quindi il risultato finale 竪 sovrapponibile. Per equalizzare il segnale sono state svolte le seguenti operazioni: 1) interpolazione numerica della curva di risposta del geofono, in modo da avere un modello matematico applicabile nei passi successivi; 2) calcolo della Trasformata di Fourier del segnale acquisito;
  • 3. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 3 di 24 3) applicazione del modello matematico numerico della curva di risposta del geofono, di cui al punto (1), nel dominio delle frequenze del segnale in elaborazione; 4) calcolo della Antitrasformata di Fourier per ricostruire il segnale corretto nel dominio del tempo. Per loperazione (1) 竪 stato utilizzato un foglio elettronico, facilmente realizzabile o con software free quali OpenOffice, LibreOffice, ecc. o commerciale quale Excel. Per le operazioni descritte ai punti (2), (3) e (4) 竪 stato utilizzato un ambiente matematico denominato Scilab. Dal sito http://www.scilab.org/ : Scilab 竪 un pacchetto di programmi gratuiti per la computazione numerica sviluppati dallo INRIA e dallo ENPC inFrancia, poi da Scilab Consortium in seno alla Fondazione Digiteo. Oggi Scilab 竪 sviluppato da Scilab Enterprises. Scilab attualmente utilizza la licenza CeCILL v2, che 竪 compatibile con la licenza GNU GPL v2 della Free Software Foundation. La versione utilizzata per le analisi 竪 la 5.5.2 sia in ambiente a 32 bit sia in quello a 64 bit (s.o. Windows 7 pro). Tale scelta 竪 stata fatta in applicazione della Mission del sistema Theremino, ossia: Ricongiungere il mondo digitale col mondo reale e concreto. Con tecnologie semplici e facilmente riproducibili. E con il massimo rispetto per lambiente e per tutti gli esseri viventi. In questo caso, mediante un acquisitore autocostruito, 竪 stato campionato un segnale relativo a un microtremore. La sua elaborazione 竪 stata finora condotta in un ambiente software specializzato, quale Geopsy. Adesso si mostra come iniziare ad operare in un ambiente software pi湛 generico, quale Scilab, rendendo ancora pi湛 evidente la possibilit di sfruttare a pieno la potenzialit dei computer, volendo anche con alternative free a Scilab quali Octave, SciPy, python(x,y), ecc. o commerciali quali Matlab, Mathematica, ecc.. In appendice viene riportato il listato delle istruzioni utilizzate in ambiente Scilab. Nel seguito si illustrano i risultati ottenuti. Nella figura seguente si riporta (1) linterpolazione della curva di risposta del geofono
  • 4. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 4 di 24 DigitalizzazionedellacurvadirispostadelgeofonoSensheSS-4.5N frequenzarispostav/(m/s)decadimentodecadimento% logXlogYXY(28,8-Y)(28,8-Y)/28,8% 0,0000,1501,0001,4127,4395,10 0,5151,1773,27315,0313,8147,88 0,5791,2803,79319,059,7933,93 0,6171,3304,14021,387,4625,87 0,6561,3704,52923,445,4018,72 0,6921,4004,92025,123,7212,90 0,7431,4305,53426,921,926,67 0,7891,4606,15228,840,000,00 0,8461,4607,01528,840,000,00 1,7161,46052,00028,840,000,00 2,3051,460201,83728,840,000,00 2,7591,460574,11628,840,000,00 2,9651,460922,57128,840,000,00 3,0001,4601000,00028,840,000,00 Interpolazionematematicapernpunticonpolinomiodigradon-1 n属XY 11,0001,413Controllointerp.012345 23,27315,031yta0a1a2a3a4a5 33,79319,055a0221,890954619,0551,00E+003,79E+001,44E+015,46E+012,07E+021,22E+03 44,14021,380a1-192,357807821,3801,00E+004,14E+001,71E+017,10E+012,94E+021,22E+03 54,52923,442a264,8499648623,4421,00E+004,53E+002,05E+019,29E+014,21E+021,91E+03 64,92025,119a3-9,30523319125,1191,00E+004,92E+002,42E+011,19E+025,86E+022,88E+03 75,53426,915a40,49615233526,9151,00E+005,53E+003,06E+011,69E+029,38E+025,19E+03 86,15228,840a5-0,00092490828,8401,00E+006,15E+003,78E+012,33E+021,43E+038,81E+03 97,01528,840 1052,00028,840 11201,83728,840 12574,11628,840 13922,57128,840 141000,00028,840 Parametri Rettada(1)a(3) Rettada(8)a(14) Matricedivalutazionedelpolinomio 1 1 2 210 ++++=n nxa...xaxaay
  • 5. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 5 di 24 ossia: (a) una retta nel piano log10-log10; (b) una interpolazione polinomiale del 5属 ordine per il tratto di raccordo; (c) un tratto di plateau con risposta costante del geofono (valore massimo 28.84 V/m/s). Per quanto riguarda i passi di calcolo (2), (3) e (4) nel seguito si illustrano i risultati dellelaborazione ottenuti in ambiente Scilab con riferimento al segnale acquisito con gain 2100 (24k Ohm). Passo di calcolo (2), la Trasformata di Fourier del segnale acquisito:
  • 6. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 6 di 24 Passo di calcolo (3), lequalizzazione della Trasformata di Fourier con linterpolazione della curva di risposta del geofono
  • 8. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 8 di 24 Passo di calcolo (4), la ricostruzione del segnale nel dominio del tempo
  • 9. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 9 di 24 Infine un confronto tra un estratto del file .saf originario e di quello ottenuto a seguito dellequalizzazione, che consente di capire la compatibilit del nuovo file con Geopsy segnale originario segnale equalizzato
  • 10. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 10 di 24 Utilizzando Geopsy si espongono i risultati dellelaborazione in termini di spettro di risposta dei singoli canali Segnale originario: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti- triggering on raw signal, n属 finestre 134 Segnale equalizzato: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti- triggering on raw signal, n属 finestre 122
  • 11. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 11 di 24 Utilizzando la funzione Spectrum Rotate si ottiene Segnale originario: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti- triggering on raw signal, campionamento 200 punti, n属 finestre 41 Segnale equalizzato: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti- triggering on raw signal, campionamento 200 punti, n属 finestre 20 Nel caso specifico si osservi che lequalizzazione delle basse frequenze poco ha aggiunto alle analisi gi svolte in termini di individuazione delle frequenze proprie dominanti del manufatto.
  • 12. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 12 di 24 Infine utilizzando la funzione H/V si ottiene Segnale originario: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti-triggering on raw signal, n属 finestre 24 Segnale equalizzato: gain 2.100 (24k ohm), lunghezza finestra da 20 a 50s, bad sample thresold 50%, anti-triggering on raw signal, n属 finestre 20 La sostanziale coincidenza dei risultati conferma quanto esposto in precedenza sulla non necessit di equalizzare un segnale per analisi H/V. Quindi adesso il processo di calcolo 竪 abbastanza completo sulla manipolazione del singolo segnale nel dominio delle frequenze. Nei successivi test si mostreranno altri aspetti di calcolo e comportamento applicati bagli edifici. Grazie per lattenzione.
  • 13. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 13 di 24 APPENDICE 1 Evoluzione dellacquisitore USB Lacquisitore Hear Our Earth ha avuto una evoluzione: 1 per ogni canale 竪 possibile selezionare un diverso valore di gain, correlato al valore di resistenza posizionata tra il geofono e lamplificatore. Le resistenze selezionabili sono 0 Ohm, 150k Ohm e 24k Ohm. La selezione avviene mediante commutatore rotativo; 2 il computer per lacquisizione 竪 un tablet Asus M81C con sistema operativo Windows 8.1, che consente lacquisizione contemporanea di almeno 6 canali a 500Hz (non avendo altri acquisitori a disposizione per un testare un maggiore carico). Per agevolare limmissione dei dati si pu嘆 utilizzare una tastiera/mouse bluetooth Ipazzport- KP-810-19BTT. Ovviamente 竪 occorso un po di lavoro per le modifiche, ma questo fa parte della volont di sperimentare.
  • 14. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 14 di 24 alcune fasi delle modifiche nota: il corto con 375 Ohm 竪 gi compreso nelle caratteristiche del geofono
  • 15. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 15 di 24 APPENDICE 2 Equalizzazione del segnale Software: Scilab ver. 5.5.2 con Windows 7 pro a 64 bit Utilizzo: sulla console di Scilab vengono visualizzati i passi di calcolo, illustrati nel seguito. In caratteri di colore nero evidenziano quanto appare a video, con quelli di colore rosso e in corsivo si riporta un commento al singolo passo di calcolo. *** FASE 1 DI 8: SELEZIONE FILE DATI *** mediante una maschera di windows 竪 possibile selezionare il file dati con estensione .saf *** INIZIO *** *** FINE *** *** FASE 2 DI 8: LETTURA FILE DATI *** vengono letti i dati dal file *** INIZIO *** *** FINE *** *** FASE 3 DI 8: CALCOLO FFT *** calcolo della Fast Fourier Trasformcon la libreria fftw3, riferimenti in http://www.fftw.org/ *** INIZIO *** *** FINE *** *** FASE 4 DI 8: VISUALIZZAZIONE FFT *** grafici dei segnali relativi ai 3 canali Verticale, Nord e Est *** INIZIO *** *** FINE *** *** FASE 5 DI 8: EQUALIZZAZIONE *** applicazione del modello numerico della curva di decadimento del geofono *** INIZIO *** *** FINE *** *** FASE 6 DI 8: RICOSTRUZIONE DEL SEGNALE *** calcolo della Antitrasformata di Fourier al segnale equalizzato nel dominio delle frequenze *** INIZIO *** *** FINE *** *** FASE 7 DI 8: SCRITTURA SEGNALE CORRETTO *** scrittura su file (con nome modificato rispetto loriginale con laggiunta in coda di _EQ.saf) *** INIZIO *** *** FINE *** *** FASE 8 DI 8: CANCELLAZIONE VARIABILI *** rilascio della memoria utilizzata *** INIZIO *** *** FINE ***
  • 16. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 16 di 24 Listato: // Equalizzazione di un segnale geofonico // Geofoni SS-4.5N della Senshe //---------------------------------- // Data: 08.05.2015 // Autore: Corrado Pecora mprintf('%s n',"*** INIZIO ELABORAZIONE ***") // Inizializzazione della memoria // Sistema Operativo: Windows 7 Pro 64 bit // Ram: 16 Gb stacksize(100000000) // 1 - SELEZIONE FILE DATI mprintf('%s n',"*** FASE 1 DI 8: SELEZIONE FILE DATI ***") // maschera di dialogo windows // selezione di un solo file SAF mprintf('%s n',"*** INIZIO ***") FileName=uigetfile(["*.saf"],"","Seleziona un file SAF",%f) mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 2 - LETTURA DATI DAL FILE mprintf('%s n',"*** FASE 2 DI 8: LETTURA FILE DATI ***") mprintf('%s n',"*** INIZIO ***") // apertura file testo in lettura fd_r=mopen(FileName,'rt') // lettura delle prime 9 righe for i=1:9 riga=mgetl(fd_r,1) end // lettura della frequenza di campionamento // (sulla stessa riga: stringa %s, stringa %s, intero %g) riga=mfscanf(fd_r,'%s %s %g') SAMP_FREQ=riga(:,3) // lettura del numero di dati // (sulla stessa riga: stringa %s, stringa %s, intero %g) riga=mfscanf(fd_r,'%s %s %g') NDAT=riga(:,3)
  • 17. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 17 di 24 // lettura di 9 righe (una per EndOfLine della precedente) for i=1:9 riga=mgetl(fd_r,1) end // lettura di terne di dati in NDAT righe relativi a velocit geofoniche segnale(1:NDAT,1:3)=0.0 for i=1:NDAT segnale(i,1:3)=mfscanf(fd_r,'%g %g %g') end // chiusura file testo mclose(fd_r) // asse dei tempi del segnale acquisito // ascisse equispaziate di dt = 1/SAMP_FREQ dt=1.0/SAMP_FREQ tempo = 0.0:dt:(NDAT-1)*dt mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 3 - FFT mprintf('%s n',"*** FASE 3 DI 8: CALCOLO FFT ***") mprintf('%s n',"*** INIZIO ***") // suddivisione del segnale segnale_V=segnale(1:NDAT,1) segnale_N=segnale(1:NDAT,2) segnale_E=segnale(1:NDAT,3) // calcolo fft_segnale_V=fft(segnale_V,-1,"symmetric") fft_segnale_N=fft(segnale_N,-1,"symmetric") fft_segnale_E=fft(segnale_E,-1,"symmetric") mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 4 - VISUALIZZAZIONE DELL'ELABORAZIONE mprintf('%s n',"*** FASE 4 DI 8: VISUALIZZAZIONE FFT ***") mprintf('%s n',"*** INIZIO ***") // essendo "segnale" un vettore di numeri reali // la fft sar coniugata simmetrica. Quindi verranno visualizzati solo
  • 18. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 18 di 24 // i primi NDAT/2 punti. // frequenze associate // VEDI http://usingscilab.blogspot.it/2009/08/fft-specifying-frequencies.html // http://www.mathworks.com/matlabcentral/answers/36430-plotting-the-frequency-spectrum dF=SAMP_FREQ/NDAT frequenza=(0:NDAT-1)*dF NrDati=size(frequenza,'*') freq_plot=frequenza(1:NrDati/2) // i 3 assi distinti VVV=fft_segnale_V(1:NrDati/2); //asse Verticale NNN=fft_segnale_N(1:NrDati/2); //asse Nord EEE=fft_segnale_E(1:NrDati/2); //asse Est // grafici (nota: frequenza' 竪 la trasposta di frequenza) scf(1) clf(1) subplot(3,1,1) plot(freq_plot',real(VVV)) xtitle("Asse Verticale - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(VVV)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(VVV)) xtitle("Valore assoluto") scf(2) clf(2) subplot(3,1,1) plot(freq_plot',real(NNN)) xtitle("Asse Nord - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(NNN)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(NNN)) xtitle("Valore assoluto")
  • 19. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 19 di 24 scf(3) clf(3) subplot(3,1,1) plot(freq_plot',real(EEE)) xtitle("Asse Est - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(EEE)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(EEE)) xtitle("Valore assoluto") mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 5 - EQUALIZZAZIONE mprintf('%s n',"*** FASE 5 DI 8: EQUALIZZAZIONE ***") mprintf('%s n',"*** INIZIO ***") // Nota: la curva di equalizzazione 竪 stata ricavata mediante // digitalizzazione del grafico della risposta in frequenza // parametri della interpolazione // 5.1) limiti in frequenza dei 3 campi: rettilineo, polinomiale, // rettilineo orizzontale // N.B.: la retta 竪 nel piano log10-log10 come la rappresentazione della // curva di risposta in frequenza del geofono freq_0=1.00 freq_1=3.79315 freq_2=6.09 // minima frequenza di equalizzazione freq_min=0.3 //max valore di risposta in frequenza ymax=28.84 // 5.2) punti della retta del campo (1) x0=log10(freq_0) x1=log10(freq_1) y0=log10(1.413) y1=log10(19.055)
  • 20. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 20 di 24 m=(y1-y0)/(x1-x0) // 5.3) parametri della interpolazione polinomiale del 5属 ordine a0=221.8909546 a1=-192.3578078 a2=64.84996486 a3=-9.305233191 a4=0.496152335 a5=-0.000924908 // 5.4) correzione della FFT per ogni singolo canale fft1_segnale_V(1:NrDati)=fft_segnale_V(1:NrDati) fft1_segnale_N(1:NrDati)=fft_segnale_N(1:NrDati) fft1_segnale_E(1:NrDati)=fft_segnale_E(1:NrDati) // ciclo di correzione dei valori i=2; //per i=1 frequenza(i)=0.0 freq_=frequenza(i) while freq_< freq_2 if freq_<= freq_min then // per questioni numeriche si equalizza la fft a partire da una certa // frequenza Coeff=1.0 else if freq_<=freq_1 then // tratto rettilineo Coeff=ymax/10.^((log10(freq_)-x0)*m+y0) else // interpolazione polinomiale Coeff=ymax/(a0+a1*freq_+a2*freq_^2+a3*freq_^3+a4*freq_^4+a5*freq_^5) end end fft1_segnale_V(i)=Coeff*fft_segnale_V(i) fft1_segnale_N(i)=Coeff*fft_segnale_N(i) fft1_segnale_E(i)=Coeff*fft_segnale_E(i) i=i+1 freq_=frequenza(i) end // i 3 assi distinti
  • 21. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 21 di 24 VVV1=fft1_segnale_V(1:NrDati/2); //asse Verticale NNN1=fft1_segnale_N(1:NrDati/2); //asse Nord EEE1=fft1_segnale_E(1:NrDati/2); //asse Est // 5.5) grafici scf(4) clf(4) subplot(3,1,1) plot(freq_plot',real(VVV1)) xtitle("Asse Verticale - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(VVV1)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(VVV1)) xtitle("Valore assoluto") scf(5) clf(5) subplot(3,1,1) plot(freq_plot',real(NNN1)) xtitle("Asse Nord - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(NNN1)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(NNN1)) xtitle("Valore assoluto") scf(6) clf(6) subplot(3,1,1) plot(freq_plot',real(EEE1)) xtitle("Asse Est - Parte Reale") subplot(3,1,2) plot(freq_plot',imag(EEE1)) xtitle("Parte Immaginaria") subplot(3,1,3) plot(freq_plot',abs(EEE1)) xtitle("Valore assoluto")
  • 22. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 22 di 24 mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 6 - RICOSTRUZIONE DEL SEGNALE mprintf('%s n',"*** FASE 6 DI 8: RICOSTRUZIONE DEL SEGNALE ***") mprintf('%s n',"*** INIZIO ***") // Calcolo della trasformata inversa di Fourier sulla base della // trasformata corretta segnale1_V=fft(fft1_segnale_V,1) //asse Verticale segnale1_N=fft(fft1_segnale_N,1) //asse Nord segnale1_E=fft(fft1_segnale_E,1) //asse Est // Grafici scf(7) clf(7) subplot(2,1,1) plot(tempo',segnale_V) xtitle("Asse Verticale - Segnale originario") subplot(2,1,2) plot(tempo',segnale1_V) xtitle("Asse Verticale - Segnale corretto") scf(8) clf(8) subplot(2,1,1) plot(tempo',segnale_N) xtitle("Asse Nord - Segnale originario") subplot(2,1,2) plot(tempo',segnale1_N) xtitle("Asse Nord - Segnale corretto") scf(9) clf(9) subplot(2,1,1) plot(tempo',segnale_E) xtitle("Asse Est - Segnale originario") subplot(2,1,2) plot(tempo',segnale1_E) xtitle("Asse Est - Segnale corretto")
  • 23. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 23 di 24 mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 7 - SCRITTURA DEL SEGNALE CORRETTO mprintf('%s n',"*** FASE 7 DI 8: SCRITTURA SEGNALE CORRETTO ***") mprintf('%s n',"*** INIZIO ***") // inserimento di"_EQ" alla fine del nome del file di input FileNameOut=strsubst(FileName,".saf","_EQ.saf") // apertura file input in lettura fd_r=mopen(FileName,'rt') // apertura del file output in scrittura fd_w=mopen(FileNameOut,'wt') // lettura e scrittura delle prime 19 righe for i=1:19 riga=mgetl(fd_r,1) mputl(riga,fd_w) end // chiusura del file di input mclose(fd_r) // scrittura di terne di valori per ogni riga for i=1:NDAT riga=sprintf('%8.3f %8.3f %8.3f ',segnale1_V(i),segnale1_N(i),segnale1_E(i)) mputl(riga,fd_w) end // chiusura del file di output mclose(fd_w) mprintf('%s n',"*** FINE ***") mprintf('%s n',"") // 8 - CANCELLAZIONE DELLE VARIABILI // per rilascio della memoria mprintf('%s n',"*** FASE 8 DI 8: CANCELLAZIONE VARIABILI ***") mprintf('%s n',"*** INIZIO ***") clear
  • 24. Corrado Pecora lisia.cp@gmail.com 08/05/2015 pagina 24 di 24 mprintf('%s n',"*** FINE ***") mprintf('%s n',"") mprintf('%s n',"*** FINE ELABORAZIONE ***") mprintf('%s n',"")