An´alise de Sinais de Ultra-som usando Decomposic¸ao …la2i.com/meel/scalassara05analise.pdf · - Um outro agradecimento a Erica, bolsista junior´ do CNPq, pelo aux´ılio na

Universidade Estadual de LondrinaCentro de Tecnologia e Urbanismo

Departamento de Engenharia Eletrica

Analise de Sinais de Ultra-som usando

Decomposicao Autorregressiva e

Rastreamento de Polos

Paulo Rogerio Scalassara

Carlos Dias MacielOrientador

Banca ExaminadoraCarlos Dias Maciel - Presidente

Ailton Akira Shinoda - UEL/LondrinaJose Carlos Pereira - USP/Sao Carlos

Dissertacao submetida ao Departamento de Engenharia Eletricada Universidade Estadual de Londrina, para preenchimento

dos pre-requisitos parciais para obtencao do tıtulo deMestre em Engenharia Eletrica.

Londrina, 10 de Marco de 2005.

ii

Dedico este trabalhoa minha noiva

Melissa.

iii

iv

Resumo

Neste trabalho, foram realizados estudos de sinais de ultra-som, tanto artificiais, gera-dos por modelos simplificados de tecidos biologicos, quanto reais, obtidos por inspecoesde ultra-som em phantoms. Esses modelos utilizados consideravam que os ecos retroes-palhados por tecidos moles, como o fıgado, eram gerados por dois tipos de espalhadores,um aleatorio, que modelava a parte difusa do tecido, e outro, regular, com distanciasdefinidas por uma propriedade do meio chamada espacamento medio dos espalhadores(MSS). Os sinais foram analisados usando decomposicao autorregressiva e rastreamentodos polos. Como o MSS pode ser caracterizado pelos picos do espectro de potencia doeco e os polos da decomposicao AR sao diretamente relacionados a esses picos, podia-seencontrar o valor do MSS sem a necessidade de estimar o PSD. Foram feitas variassimulacoes utilizando o metodo de Monte Carlo para encontrar os efeitos de variacoesdos parametros: jitter, que mede a regularidade do meio, amplitude da parte difusa eordem da decomposicao. Percebeu-se que o aumento do jitter ou da amplitude provocaum aumento do erro de estimativa do MSS, sendo o primeiro o parametro mais sensıvel.Sobre a ordem da decomposicao, nao foram obtidos resultados conclusivos. Os resulta-dos com sinais de phantoms foram muito bons, pois os erros eram compatıveis ou atemenores que os apresentados nas simulacoes. A grande vantagem do metodo usado ea sua simplicidade, pois nao e necessaria a estimativa do PSD e os resultados encon-trados sao similares aos da literatura. Ele tambem traz uma contribuicao aos estudosde ecos de ultra-som retroespalhados por meios pela analise das amplitudes relativasdos primeiros picos do PSD e suas larguras de banda quando se varia a regularidadedo meio.

v

vi

Abstract

In this work, two kinds of ultrasound signals were studied: artificial ones, generatedby simplified biological tissue models, and real ones, obtained by ultrasound inspectionof phantoms. These models considered the fact that the echoes backscattered by softtissues, like liver, were generated by two kinds of scatterers: first, a random one,modelling the diffuse part of the tissue, and second, a regular one, with positionsdefined by a tissue property called mean scatterer spacing (MSS). The signals wereanalyzed using an autoregressive decomposition and pole tracking. Since the MSSmay be characterized by the peaks of the echo power spectral density and the polesof the AR decomposition are directly related to these peaks, it was possible to findthe MSS value without estimating the PSD. Several simulations using the Monte Carlomethod were performed seeking for the effects of the following parameters variations:jitter, that measures the medium regularity, amplitude of the diffuse components anddecomposition order. It could be noticed that a raise in the jitter or the amplituderesulted in an increase of the MSS estimated error, being the jitter the most sensitiveparameter. The results for the order of the decomposition were not conclusive. Theresults with phantom signals were very good, because the errors were compatible oreven smaller than those obtained in the simulations. The best advantage of this methodis its simplicity, because it not necessary to estimate the PSD and the obtained resultsare similar to those presented by the literature. It also brings a contribution to thestudy of ultrasound backscattered echoes by the analysis of the relative amplitudes ofthe PSD first echoes and their bandwidths when the medium regularity is varied.

vii

viii

Agradecimentos

- Eu gostaria de agradecer primeiramente a Deus, por ter me guiado e iluminadoneste caminho varias vezes arduo. Tambem por ter me capacitado a nao perderas esperancas nos momentos de erros.

- Em seguida, aos meus pais por terem dado suporte neste perıodo e nos anteriorese, tambem, ajudado a levantar meu animo quando este estava para baixo.

- Um agradecimento especial para a minha noiva Melissa, por ter me incentivado,mesmo quando eu deixava de dar atencao a ela para me dedicar ao trabalho.

- Agradeco ao professor Carlos Maciel por ter me orientado e gasto horas discu-tindo os problemas encontrados, alem de ter feito de tudo para que este trabalhoseguisse por um caminho interessante para a Ciencia.

- Quero agradecer aos meus colegas pelo apoio e ajuda quando eu precisava, etambem para meus amigos da USP por terem me ajudado durante o tempo quepassei em Sao Carlos.

- Um outro agradecimento a Erica, bolsista junior do CNPq, pelo auxılio na fa-bricacao dos phantoms.

- Por fim, gostaria de agradecer aos professores que me acompanharam nesteperıodo, inclusive aqueles que, mesmo sem dar aulas, conseguiram passar umpouco de seus conhecimentos. Tambem, a UEL, pela estrutura e servicos, e aoCNPq pelo financiamento.

ix

x

Sumario

Resumo vi

Abstract viii

Agradecimentos x

Sumario xii

Lista de Figuras xix

Lista de Tabelas xxi

1 Introducao 11.1 Objetivos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.2 Justificativa . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.3 Introducao ao Ultra-som . . . . . . . . . . . . . . . . . . . . . . . . . . 2

2 Revisao Bibliografica 152.1 Modelamento de Tecidos Biologicos . . . . . . . . . . . . . . . . . . . . 152.2 Espacamento Medio dos Espalhadores . . . . . . . . . . . . . . . . . . . 19

3 Teoria 253.1 Processos Estocasticos . . . . . . . . . . . . . . . . . . . . . . . . . . . 253.2 Modelos Autorregressivos . . . . . . . . . . . . . . . . . . . . . . . . . . 303.3 Filtragem Adaptativa . . . . . . . . . . . . . . . . . . . . . . . . . . . . 403.4 Rastreamento de Polos . . . . . . . . . . . . . . . . . . . . . . . . . . . 513.5 Simulacao de Monte Carlo . . . . . . . . . . . . . . . . . . . . . . . . . 55

4 Materiais e Metodos 594.1 Phantoms . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 594.2 Equipamentos e Software de Aquisicao . . . . . . . . . . . . . . . . . . 604.3 Softwares de Simulacao e Analise . . . . . . . . . . . . . . . . . . . . . 674.4 Aspectos Computacionais . . . . . . . . . . . . . . . . . . . . . . . . . 70

5 Resultados 815.1 Resultados com Sinais Artificiais . . . . . . . . . . . . . . . . . . . . . 815.2 Simulacoes de Monte Carlo . . . . . . . . . . . . . . . . . . . . . . . . . 875.3 Resultados de Sinais Coletados de Phantoms . . . . . . . . . . . . . . . 104

xi

SUMARIO

6 Discussoes e Conclusao 1236.1 Discussao . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1236.2 Conclusao . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1306.3 Proximos Passos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 132

Referencias Bibliograficas 133

xii

Lista de Figuras

1.1 Ilustracao do espalhamento causado pela incidencia de uma onda deultra-som nos tres tipos de distribuicoes de espalhadores. As setas emlinha cheia correspondem a propagacao sem atenuacao e as linhas ponti-lhadas, com atenuacao. (a) Distribuicao tenue. (b) Distribuicao menostenue. (c) Distribuicao densa. . . . . . . . . . . . . . . . . . . . . . . . 6

1.2 Ilustracao de um trem de pulsos de ultra-som tıpico usado nas tecnicasde pulso-eco. Apresentam-se a duracao do pulso, τ , e o perıodo derepeticao do pulso, PRP [7]. . . . . . . . . . . . . . . . . . . . . . . . . 7

1.3 Ilustracao das resolucoes axial e lateral de dois sistemas. As imagensnomeadas 1 sao de um sistema com resolucoes suficientes e as 2 de umsistema com resolucoes insuficientes. (a) Resolucao axial. (b) Resolucaolateral. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

1.4 Exemplo do processo de obtencao de ecos de ultra-som de tecidos. (a)Interacao da onda de ultra-som com as interfaces do tecido mostrando asdiferencas nas reflexoes. (b) Sinal de eco recebido em funcao da distanciadas interfaces, desconsiderando o espalhamento do tecido. . . . . . . . . 9

1.5 Exemplo da Figura 1.4 considerando-se tambem o espalhamento dotecido. As interfaces sao identificadas pelos picos e os valores inter-mediarios sao referentes a constituicao do meio. . . . . . . . . . . . . . 10

1.6 Diagrama simplificado de um equipamento modo A. (a) Diagrama deblocos. (b) Exemplos dos principais sinais. . . . . . . . . . . . . . . . . 11

1.7 Exemplo de formacao de imagem modo B, utilizando-se os ecos paramodificar o brilho dos feixes de eletrons do monitor CRT. . . . . . . . . 12

1.8 Diagrama de blocos simplificado de um equipamento modo B. . . . . . 13

2.1 Ilustracao do modelo unidimensional apresentado em [21] considerandoos tecidos biologicos como uma distribuicao aleatoria de espalhadores(scatterers) com forma e forca de espalhamento aleatorios. . . . . . . . 17

2.2 (a) Esquema da constituicao do tecido de fıgado humano. (b) Corte deum fıgado humano mostrando as estruturas que o compoem (Figurasretiradas de [38]). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18

2.3 Ilustracao do espacamento medio dos espalhadores, considerando o mo-delo unidimensional apresentado anteriormente. . . . . . . . . . . . . . 19

2.4 Exemplos da resposta impulsiva do tecido usando o modelo, junto comseus espectros de potencia, para jitter de (a) 5% e (b) 50%. A variacaoda amplitude dos ecos, nestes exemplos, e imperceptıvel visualmente. . 21

xiii

LISTA DE FIGURAS

2.5 (a) Espectro de potencia da resposta do transdutor. (b) Espectro depotencia do eco mostrando o deslocamento provocado pela frequenciade ressonancia e largura de banda do transdutor. . . . . . . . . . . . . 22

2.6 (a) Espectro de potencia do exemplo de sinal de eco apresentado na fi-gura anterior. (b) Espectro de potencia do sinal ao quadrado, mostrandoo primeiro harmonico com maior amplitude, correspondendo a fMSS, ea diminuicao do efeito da resposta do transdutor. . . . . . . . . . . . . 23

3.1 Ilustracao da formacao de um processo estocastico, ~x(n), considerando-se o espaco amostral Ω, cujos eventos ωi geram as variaveis aleatoriasxk(i). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27

3.2 Diagrama de um modelo ARMA para um processo estocastico, consi-derando a entrada x[n], tipicamente ruıdo branco, e a saıda y[n]. Parailustracao, considera-se p = q. . . . . . . . . . . . . . . . . . . . . . . . 31

3.3 Diagrama de um modelo media movel (MA) para um processo estocastico,considerando a entrada x[n] e a saıda y[n]. . . . . . . . . . . . . . . . . 33

3.4 Diagrama de um modelo autorregressivo (AR) para um processo es-tocastico, considerando a entrada x[n] e a saıda y[n]. . . . . . . . . . . 33

3.5 Efeito da variacao do raio dos polos complexos conjugados de um filtroAR de segunda ordem, percebendo-se o aumento do pico e a diminuicaoda largura de banda proporcionalmente a aproximacao do raio ao cırculounitario. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35

3.6 Efeito da variacao do raio dos zeros complexos conjugados de um fil-tro MA de segunda ordem, percebendo-se o aumento da profundidadedo vale e a diminuicao da largura de banda proporcionalmente a apro-ximacao do raio ao cırculo unitario. . . . . . . . . . . . . . . . . . . . . 35

3.7 Comparacao do PSD do modelamento de um processo AR de segundaordem com banda larga (raio dos polos complexos conjugados iguais a0,7071) por modelos MA de ordem 2, 5 e 10. . . . . . . . . . . . . . . . 37

3.8 Comparacao do PSD do modelamento de um processo AR de segundaordem com banda estreita (raio dos polos complexos conjugados iguaisa 0,9592) por modelos MA de ordem 2, 5 e 10. . . . . . . . . . . . . . . 38

3.9 Aproximacao de um processo ARMA(2,2), que corresponde a um AR(2)mais ruıdo branco, por modelos AR de segunda e quinta ordens. Oszeros daquele processo sao 0, 5 e0,6 π e 0, 5 e−0,6 π, com modulos longe docırculo unitario, ou seja, a variancia do ruıdo e baixa. . . . . . . . . . . 39

3.10 Aproximacao de um processo ARMA(2,2), que corresponde a um AR(2)mais ruıdo branco, por modelos AR de segunda e quinta ordens. Oszeros daquele processo sao 0, 9 e0,6 π e 0, 9 e−0,6 π, com modulos perto docırculo unitario, ou seja, a variancia do ruıdo e alta. . . . . . . . . . . . 39

3.11 Diagrama de blocos geral dos filtros Wiener, sendo x[n] a entrada, y[n]a saıda, d[n] a resposta desejada, e[n] o erro de estimativa e w os coefi-cientes do filtro. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40

3.12 Diagrama de blocos de um preditor linear para frente, sendo x[n −1], . . . , x[n−M ] as entradas e x[n] a estimativa de x[n]. . . . . . . . . . 42

xiv

LISTA DE FIGURAS

3.13 Diagrama de blocos de um preditor linear para tras, sendo x[n], . . . , x[n−M + 1] as entradas e x[n−M ] a estimativa de x[n−M ]. . . . . . . . . 43

3.14 Valores dos coeficientes para o modelo LMS aplicado a um processoautorregressivo de segunda ordem. O grafico de cima apresenta os re-sultados para passo dividido por 50 e, o de baixo, para passo divididopor 100. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47

3.15 Diagrama de blocos do metodo de prony, mostrando a forma de se en-contrar o erro e[n], sendo b[n] a estimativa de b[n]. . . . . . . . . . . . . 48

3.16 Polos de H(z) do exemplo apresentados no plano z, sendo 0, 3626; −0, 8325+0, 4425i; −0, 8325− 0, 4425i; 0, 5012 + 0, 6077i e 0, 5012− 0, 6077i. . . . 53

3.17 Densidade espectral de potencia de y[n] mostrando dois picos, o menorcom frequencia central normalizada 0,1402 Hz, o maior com frequencia0,4222 Hz. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53

3.18 Histograma das distancias dos pontos aleatorios ate a origem para o casode n = 500. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56

3.19 Erro de estimativa para 100 simulacoes com diferentes valores de n,comecando em 500 ate 50000. Percebe-se o erro diminui com o aumentodo numero de pontos. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57

4.1 Ilustracao da constituicao dos phantoms usados como meios de testepara inspecao de ultra-som. . . . . . . . . . . . . . . . . . . . . . . . . 60

4.2 Foto de tres phantoms usados para se coletar sinais. Pode-se ver asmarcas de tres interfaces no phantom da esquerda e quatro no da direita,no do meio nao e possıvel se observar. . . . . . . . . . . . . . . . . . . . 60

4.3 Os dois transdutores usados nas coletas de sinais de ultra-som. O dadireita e focalizado com diametro de 10 mm e frequencia de operacaode 8 MHz, e o da esquerda e nao focalizado com diametro de 15 mm efrequencia de 3, 5 MHz. . . . . . . . . . . . . . . . . . . . . . . . . . . 61

4.4 Ilustracao da montagem experimental usada nas coletas de sinais deultra-som de meios de teste. . . . . . . . . . . . . . . . . . . . . . . . . 61

4.5 Foto do equipamento usado nos testes, montado no Laboratorio de Ins-trumentacao Biomedica da UEL. . . . . . . . . . . . . . . . . . . . . . 62

4.6 Ilustracao da disposicao do campo de ultra-som criado por um transdu-tor circular com impedancia casada com a do meio. . . . . . . . . . . . 63

4.7 Ilustracao dos campos proximo e distante para um transdutor focalizadono qual se utiliza lente, apresentando-se tambem a regiao de foco. . . . 64

4.8 Interface do programa de controle e aquisicao de dados, SADUS, atual-mente na versao 5.0. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65

4.9 Modelo geral utilizado pelo software para salvar os dados em arquivos.Apresenta-se tambem o cabecalho usado para identificacao dos dados [42]. 66

4.10 Codificacao usada no cabecalho dos arquivos para identificacao dos da-dos. A tensao e dada em volts, a frequencia em hertz. . . . . . . . . . . 66

4.11 Ilustracao dos resultados do algoritmo gerador de ecos. (a) Pulso deultra-som. (b) PSD do pulso. (c) Resposta impulsiva do meio. (d) Ecosimulado. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68

xv

LISTA DE FIGURAS

4.12 (a) Envoltoria de uma janela de 512 pontos do eco apresentado na Fi-gura 4.11, com Ad = 15% e V = 10%. (b) PSD desse sinal janeladomostrando os picos das tres primeiras harmonicas. . . . . . . . . . . . . 70

4.13 Graficos com os modulos dos polos da decomposicao AR em funcaoda frequencia (ate 2 MHz) utilizando os dois metodos propostos. (a)Metodo da autocorrelacao. (b) Algoritmo LMS. . . . . . . . . . . . . . 71

5.1 Sinal de envoltoria com media nula usado na analise do caso de teste 1(em cima) e espectro de potencia do sinal (embaixo). . . . . . . . . . . 82




5.5 Modulo dos polos da decomposicao AR do sinal usado no caso 1 emfuncao da frequencia (somente a parte positiva do espectro). . . . . . . 85




5.9 Variacoes dos valores de erro para os quatro casos de teste. . . . . . . . 88

5.10 Histograma das ocorrencias dos valores de MSS encontrados para jitterigual a 5% para a primeira harmonica. (a) Ad = 15%. (b) Ad = 30%. . 88

5.11 Histograma das ocorrencias dos valores de MSS encontrados para jitterigual a 15% para a primeira harmonica. (a) Ad = 15%. (b) Ad = 30%. 89



5.14 Resultados consolidados para as simulacoes de variacao de jitter mos-trando os erros de MSS para as tres harmonicas. (a) Ad = 15%. (b) Ad= 30%. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90

5.15 Resultados consolidados para as simulacoes com variacao de quantidadede ruıdo, mostram-se os erros de MSS para as tres harmonicas. (a) Jitter= 10%. (b) Jitter = 30%. . . . . . . . . . . . . . . . . . . . . . . . . . 91

5.16 Histograma das ocorrencias dos valores de MSS encontrados para Adigual a 10% para a primeira harmonica. (a) Jitter = 10%. (b) Jitter =30%. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92


xvi

LISTA DE FIGURAS



5.20 Percentual de acertos medios do MSS para a primeira harmonica em 512simulacoes variando-se o jitter e Ad. (a) Margem de 5%. (b) Margemde 25%. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94

5.21 Percentual de acertos medios do MSS para a segunda harmonica em 512simulacoes variando-se o jitter e Ad. (a) Margem de 5%. (b) Margemde 25%. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95

5.22 Percentual de acertos medios do MSS para a terceira harmonica em 512simulacoes variando-se o jitter e Ad. (a) Margem de 5%. (b) Margemde 25%. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96

5.23 Valores obtidos para os parametros jitter e Ad utilizando a regra apre-sentada. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96

5.24 Erro medio do MSS para a primeira harmonica em 512 simulacoes com1024 combinacoes de valores de jitter e Ad. . . . . . . . . . . . . . . . . 97

5.25 Erro medio do MSS para a segunda harmonica em 512 simulacoes com1024 combinacoes de valores de jitter e Ad. . . . . . . . . . . . . . . . . 97

5.26 Erro medio do MSS para a terceira harmonica em 512 simulacoes com1024 combinacoes de valores de jitter e Ad. . . . . . . . . . . . . . . . . 98

5.27 Erro medio do MSS para a primeira harmonica com 1024 combinacoesdos parametros considerando-se somente polos com frequencias menoresque 2 MHz. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98

5.28 Erro medio do MSS para a segunda harmonica com 1024 combinacoesdos parametros considerando-se somente polos com frequencias menoresque 2 MHz. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99

5.29 Erro medio do MSS para a terceira harmonica com 1024 combinacoesdos parametros considerando-se somente polos com frequencias menoresque 2 MHz. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99

5.30 Distancia acumulada dos modulos dos tres maiores polos das simulacoesao cırculo unitario em funcao das variacoes de jitter e Ad. . . . . . . . 100

5.31 Erro medio do MSS de 1,25 mm, equivalente a 40 amostras, para 128simulacoes, variando-se a ordem da decomposicao AR de 5 a 100, de 5em 5 unidades. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101


5.33 Erro medio do MSS de 1,5 mm, equivalente a 49 amostras, para 128simulacoes, variando-se a ordem da decomposicao AR de 32 a 73, de 1em 1 unidade. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102


xvii

LISTA DE FIGURAS

5.35 Erro medio do MSS de 2,5 mm, equivalente a 80 amostras, para 128simulacoes, variando-se a ordem da decomposicao AR de 25 a 225, de25 em 25 unidades. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103



5.38 Sinais coletados do phantom no. 1 com transdutor de 3,5 MHz, mostram-se o sinal de referencia (em cima) e o da amostra (embaixo). . . . . . . 107

5.39 Sinal do phantom no.1, N = 512, taxa de amostragem de 12,5 MHz (emcima). Densidade espectral de potencia da envoltoria do sinal (embaixo). 107

5.40 Polos da decomposicao AR para o sinal do phantom no.1, apresentam-seos polos de fase positiva com frequencias menores do que 1 MHz. . . . 108


5.42 Polos da decomposicao AR do sinal do phantom no.2, polos com frequenciasmenores do que 1 MHz. (a) Ordem igual a 57. (b) Ordem igual a 68. . 109


5.44 Polos da decomposicao AR para o sinal do phantom no.3, polos comfrequencias menores do que 1 MHz. (a) Ordem igual a 63. (b) Ordemigual a 96. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111


5.46 Polos da decomposicao AR para o sinal do phantom no.4, apresentam-seos polos de fase positiva com frequencias menores do que 1 MHz. (a)Ordem igual a 64. (b) Ordem igual a 78. . . . . . . . . . . . . . . . . . 112






5.52 Polos da decomposicao AR para o sinal do phantom no.7, apresentam-seos polos de fase positiva com frequencias menores do que 1 MHz, ordemigual a 88. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 117


xviii

5.54 Polos da decomposicao AR para o sinal do phantom no.8, polos comfrequencias menores do que 1 MHz. (a) Ordem igual a 107. (b) Ordemigual a 103. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119



5.57 Sinal do phantom no.10, N = 512, taxa de amostragem de 12,5 MHz(em cima). Densidade espectral de potencia da envoltoria do sinal (em-baixo). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 121


xix

LISTA DE FIGURAS

xx

Lista de Tabelas

1.1 Velocidade de propagacao do som, impedancia caracterıstica e coefici-ente de atenuacao para varios meios conhecidos [11]. . . . . . . . . . . . 4

4.1 Parametros de aquisicao de dados do programa SADUS. . . . . . . . . 64

5.1 Parametros mantidos constantes nas analises dos quatro casos de teste. 825.2 Parametros variados nas analises dos quatro casos de teste. . . . . . . . 825.3 Resultados encontrados para os quatro casos de teste. . . . . . . . . . . 875.4 Resultados medios encontrados em 100 simulacoes para os quatro casos

de teste. Apresentam-se a media (µ) e desvio padrao (σ) dos resultados. 87

xxi

LISTA DE TABELAS

xxii

Capıtulo 1

Introducao

Neste trabalho, realiza-se um estudo sobre tecidos biologicos utilizando inspecao porultra-som em meios artificiais de teste, chamados phantoms, os quais simulam as carac-terısticas de tecidos reais. O enfoque principal do trabalho e o processamento digitaldos ecos capturados das inspecoes dos meios e tambem em sinais criados por modelosmatematicos utilizando tecnicas de parametrizacao e rastreamento.

Apresentam-se os equipamentos e softwares utilizados nos experimentos comultra-som junto com a descricao dos meios e sinais envolvidos, alem de uma revisaoda teoria envolvida. A abordagem utilizada e a decomposicao autorregressiva des-ses sinais, busca-se avaliar a localizacao dos polos do sistema e monitorar suas movi-mentacoes quando parametros dos meios, como o espacamento medio dos espalhadoresou a quantidade de ruıdo presente no sinal, sao modificados.

1.1 Objetivos

O objetivo deste trabalho e realizar estudos de sinais de ultra-som simulados e coletadosde meios de teste utilizando metodos de processamento de sinais nao extensivamenteexplorados na literatura. Pela analise dos sinais, objetiva-se encontrar um parametrodo meio chamado espacamento medio dos espalhadores sem a necessidade de se estimara densidade espectral de potencia do eco.

Assim, deseja-se encontrar a eficiencia da tecnica em estimar esse parametrotanto para sinais simulados quanto para reais. E, a partir disso, observar os efeitos davariacao do nıvel de ruıdo do sinal e da regularidade do meio nos resultados obtidos.

Por fim, tem-se por objetivo encontrar uma relacao entre a diminuicao das am-plitudes das harmonicas do espectro de potencia (obtidas pelos polos da decomposicaoautorregressiva do sinal) e o aumento de suas larguras de banda com o grau de regu-laridade e/ou quantidade de ruıdo do eco.

1.2 Justificativa

Faz-se um estudo por analise de ultra-som, pois esse metodo apresenta varias vantagensem relacao a outros, como o Raio-X. Pode-se citar a possibilidade de produzir imagensem tempo real e nao ser uma forma de radiacao ionizante. Utiliza-se a inspecao modo

1

Capıtulo 1 - Introducao

A, que fornece informacoes qualitativas das estruturas do meio, as quais podem serusadas para correlacao com patologias pela comparacao entre os resultados de tecidosnormais e doentes em conjunto com analises biologicas. Essa area de pesquisa e muitopromissora, devido a quantidade de estudos e publicacoes recentes.

Os estudos sobre propriedades dos tecidos como a regularidade das estruturassao importantes, pois, conforme apresentado em trabalhos recentes, a variacao dessesparametros pode indicar problemas ou mesmo doencas, como e o caso da cirrose emtecidos de fıgado.

1.3 Introducao ao Ultra-som

A aplicacao do ultra-som para diagnosticos medicos era conhecida desde a decada de1940, quando se comecou a explorar suas capacidades utilizando a tecnologia do sonardesenvolvida durante a Segunda Guerra Mundial. Nos Estados Unidos, um grupoformado pelo medico John Wild e o engenheiro John Reid construiu um instrumentode pulso-eco, demonstrando que o ultra-som podia ser usado para detectar tumores emvarios tipos de tecidos. Outras frentes de estudo, como a formada por Willian Fly eFloyd Dunn na Universidade de Illinois, mostraram que o ultra-som de alta intensidadepodia ser usado para fins terapeuticos, quando conseguiram modificar estruturas docerebro de pacientes no tratamento do mal de Parkinson [7].

Embora tenha sido uma ferramenta utilizada em muitos diagnosticos medicosdesde aquela epoca, so foi totalmente aceita a partir da decada de 1970, quando secomecou a usar tons de cinza em imagens geradas por ultra-som para mostrar o graude reflexao de um meio. Desde entao, suas aplicacoes cresceram bastante, principal-mente devido a algumas caracterısticas vantajosas, como: ser uma forma de radiacaonao-ionizante, considerada segura pelo conhecimento cientıfico atual; produzir imagensem tempo real; ter uma resolucao de alcance milimetrico para as frequencias usadasatualmente; e ser utilizado para se estudar o fluxo de lıquidos pelo princıpio Doppler[7].

O princıpio fısico para a geracao de ondas de ultra-som e o efeito piezeletrico, oqual e uma propriedade de alguns tipos especiais de cristais, que reagem a aplicacaode uma pressao mecanica produzindo uma polarizacao eletrica proporcional, sendoque o efeito oposto tambem e observado. A piezeletricidade foi descoberta por Pierre eJacques Curie em 1880, sendo piezo uma palavra derivada do grego piezin, que significapressionar [38]. Outros materiais com essa propriedade, como as ceramicas (materiaisferroeletricos), podem ser criados por processos industriais, sendo o PZT (zirconatotitanato de chumbo) o mais comum, devido ao seu forte efeito piezeletrico. Utilizam-seas ceramicas na construcao de transdutores, os quais sao dispositivos de convercao deuma forma de energia em outra, sendo que, nas aplicacoes de ultra-som, usam-se ostransdutores piezeletricos.

As propriedades acusticas sao importantes para uma boa utilizacao do ultra-somno estudo de materiais. Entre algumas, pode-se citar: a velocidade do som no meio,a impedancia acustica, a compressibilidade e os coeficientes de atenuacao, retroespa-lhamento e reflexao. Essas caracterısticas sao relevantes porque podem ser obtidaspor tecnicas nao-invasivas de ultra-som, sem a necessidade da criacao de imagens, pois

2


possuem mais informacoes sobre o meio inspecionado [17, 22]. Em tecidos biologicos,a analise desses parametros quantitativos e suas variacoes pode indicar mudancas nostecidos, aumentando, de forma significativa, a utilizacao do ultra-som na determinacaode doencas [3, 18, 22].

Algumas das caracterısticas vantajosas do uso do ultra-som como ferramenta dediagnostico ja foram citadas, como a nao-ionizacao do meio em que se propaga, por seruma onda mecanica. Mas, alem disso, outro detalhe importante e nao haver translacaode partıculas individuais, pois o ultra-som e gerado por uma deformacao em umaceramica, assim, a energia mecanica e transferida pela oscilacao das partıculas.

Se as ondas de ultra-som possuem baixa intensidade, estas passam pelos tecidosvivos sem alterar suas funcoes, mas as de alta energia podem produzir aquecimentoe cavitacao, o que acarretaria danos celulares. O aquecimento e causado pela energiaatenuada e espalhada pela passagem da onda no tecido; e a cavitacao, que ocorre adensidades altas de potencia, e o processo de aparecimento de bolhas de gas devido aoalto gradiente de pressao local [38].

Alguns conceitos basicos de propagacao acustica devem ser apresentados, como ocomprimento de onda, λ, com unidade em metros (m), definido como a distancia entreos pontos de maxima (ou mınima) amplitude da oscilacao das partıculas do meio ondea onda se propaga. A frequencia da onda, f , expressa em Hertz (Hz); o perıodo, T ,em segundos (s), que e o inverso da frequencia. Essas grandezas sao relacionadas pelaEquacao 1.1, sendo c a velocidade de propagacao, em metros por segundo (m/s).

c = λf (1.1)

A velocidade do som, que e diretamente proporcional a rigidez do meio, dependesomente da temperatura e desse meio em que se propaga (se este for isotropico), sendodeterminada pelo atraso entre o movimento de partıculas adjacentes, o qual dependeda compressibilidade, G (em m2/N), e da densidade do meio, ρ (em kg/m3). Essarelacao e apresentada na Equacao 1.2, [7]. Para aplicacoes medicas, a velocidade emtecidos moles, como os tecidos do fıgado, pode ser considerada constante e igual a 1540m/s [11, 38].

c =1√ρG

(1.2)

A impedancia acustica caracterıstica de um meio e definida como a razao entrea pressao e a velocidade de propagacao do meio, de forma analoga a tensao e correnteem um circuito eletrico, portanto, com o mesmo significado fısico. Para fluidos, aimpedancia acustica, Z, com unidade kg/m2s ou Rayl, e apresentada pela Equacao1.3, [7, 9].

Z = ρc (1.3)

O comprimento de onda do ultra-som nas frequencias usadas em aplicacoes medicastem dimensoes comparaveis a de muitas estruturas dos tecidos, assim, as interacoes en-tre estas acontecem de maneiras complexas [9]. Quando a onda penetra em um tecido,a sua energia e absorvida por ele, sendo convertida em calor e espalhada para outras

3


direcoes sempre que encontra uma descontinuidade nas caracterısticas acusticas. Por-tanto, ao atravessar um tecido uniforme, a intensidade da onda decresce ao longo dotrajeto, sendo essa diminuicao chamada atenuacao, a qual apresenta uma dependenciada frequencia aproximadamente linear para tecidos moles na faixa de 1 a 50 MHz [7].

Para se quantificar essa caracterıstica, utiliza-se o coeficiente de atenuacao, oqual e uma medida usada para descrever a habilidade de um meio em atenuar ondas deultra-som. Sua unidade e o dB/m. O valor para tecidos moles e de aproximadamente0, 3 dB/cmMHz [7]. Para uma dada frequencia, pode-se calcular o coeficiente deatenuacao, α, utilizando-se a amplitude incidente (Ai) e emergente (Ae) conforme aEquacao 1.4, sendo z a distancia percorrida pela onda. A Tabela 1.1 apresenta avelocidade de propagacao do som, impedancia caracterıstica e coeficiente de atenuacaopara varios meios conhecidos.

α =20 log(Ae/Ai)

z(1.4)

Tabela 1.1: Velocidade de propagacao do som, impedancia caracterıstica e coeficientede atenuacao para varios meios conhecidos [11].

Meio Velocidade de Impedancia Coeficiente dePropagacao Acustica Atenuacao

(m/s) (MRayl) a 1 MHz (dB/cm)Ar 330 0,0004 10

Pulmao 650 - 1160 0,26 - 0,46 40Agua 1480 1,52 0,002

Tecidos Moles 1540 1,35 - 1,68 0,3 - 1,5Sangue 1550 1,61 0,18

Osso (Cranio) 3500 6 10Alumınio 6400 17 0,02

Quando uma onda de ultra-som que esteja se propagando por um meio encontraum segundo com incidencia normal (perpendicular), parte de sua energia e refletida eoutra continua sem desvio. A onda refletida e chamada de eco, e depende das carac-terısticas dos dois meios. Para uma onda plana, o coeficiente de reflexao R e definidopela Equacao 1.5, sendo Z1 e Z2 as impedancias acusticas desses meios respectivamente[20].

R =

(Z2 − Z1

Z2 + Z1

)2

(1.5)

Quando as impedancias dos meios sao iguais, o coeficiente de reflexao e zero,portanto toda a energia e transmitida. Ao contrario, quando as impedancias sao muitodiferentes, o coeficiente e bem proximo da unidade, entao quase toda a energia e re-fletida. A diferenca de impedancias tem implicacoes importantes no uso do ultra-som,pois sabe-se que uma estrutura abaixo de tecidos com muito ar, como o pulmao, nao edetectada pela inspecao de ultra-som, devido as baixas impedancias. Isso, pois a ondae quase totalmente refletida pelo tecido com ar, nao chegando ate a outra estrutura.

4


Se a onda de ultra-som nao incidir perpendicularmente a interface, ocorre di-vergencia, sendo o angulo de incidencia e de transmissao relacionados pela chamadaLei de Snell [20]. Quando a interface e rugosa, o que e comum na pratica, a reflexao echamada difusa, em oposicao a reflexao especular da interface plana.

A atenuacao nao e so causada pela absorcao da onda, mas tambem, pelo espalha-mento, pois, a energia da onda pode se espalhar quando entra em um meio, conforme forencontrando pequenas falhas na homogeneidade, ocupando, assim, uma area maior aoavancar. O processo de espalhamento de ultra-som em tecidos biologicos e um fenomenocomplicado, mas como este e de muita importancia para a caracterizacao desses meios,varios estudos foram feitos nessa area [38]. Estes tratam os tecidos de duas maneirasdiferentes: a primeira como uma distribuicao de espalhadores (scatterers) discretos e,segundo, como um contınuo com variacoes de densidade e compressibilidade.

A forma utilizada neste trabalho e a primeira, descrita extensivamente em [16].Neste, os tecidos biologicos sao considerados como meios aleatorios, sendo que as ondasque neles se propagam apresentam variacoes randomicas em amplitude e fase, devendoser tratadas em termos estatısticos e probabilısticos. A abordagem para analise doespalhamento e dividida em duas etapas: primeira, considera-se o espalhamento eabsorcao de um unico espalhador, e segundo, considera-se as caracterısticas de ondaquando muitos espalhadores sao distribuıdos aleatoriamente.

Para a analise de um unico espalhador, faz-se uso de uma de suas propriedadesfısicas, a secao reta de espalhamento, a qual relaciona as densidades de fluxo de potenciadas ondas incidente e espalhada em uma dada direcao. Essa caracterıstica depende,de maneira geral, da capacidade de espalhamento da partıcula e das direcoes de in-cidencia e espalhamento. Pode-se usar tambem a secao reta de retroespalhamento, coma consideracao da onda refletida na direcao da incidente, e a secao reta de absorcao,que mede a capacidade de absorcao do espalhador (tambem dependente da direcao deincidencia).

Para a segunda etapa da analise, quando a onda se propaga em meios contendovarias partıculas, considera-se dois casos: distribuicoes tenues e densas. No primeiro,pode-se usar a aproximacao de espalhamento singular, pois, assume-se que a onda inci-dente do transmissor alcanca o receptor apos encontrar poucas partıculas, desprezando-se multiplos espalhamentos. Nessa aproximacao, somam-se as potencias espalhadas porcada partıcula considerada.

Em distribuicoes com maiores densidades (mas ainda tenues), deve-se considerara atenuacao por espalhamento e absorcao pelo caminho. Essa aproximacao e chamadade multiplo espalhamento de primeira ordem. Para distribuicoes densas, utilizam-seaproximacoes de espalhamento multiplo ou de difusao. A Figura 1.1 ilustra o espa-lhamento para esses tres tipos de distribuicao apresentados, respectivamente os items(a), (b) e (c). As setas em linha cheia correspondem a propagacao sem atenuacao e aslinhas pontilhadas, com atenuacao.

Uma consideracao importante e a questao da coerencia do campo. Como aspartıculas de um meio sao distribuıdas aleatoriamente, o campo espalhado nao e cons-tante e sua amplitude e fase devem flutuar de maneira randomica, assim, a potenciatambem flutua dessa forma. Pode-se escrever o campo como uma soma de duas compo-nentes: uma media, chamada coerente, e outra aleatoria, chamada incoerente. Assim,em analises de espalhamento, pode-se considerar que as ondas espalhadas por distri-

5


Figura 1.1: Ilustracao do espalhamento causado pela incidencia de uma onda de ultra-som nos tres tipos de distribuicoes de espalhadores. As setas em linha cheia corres-pondem a propagacao sem atenuacao e as linhas pontilhadas, com atenuacao. (a)Distribuicao tenue. (b) Distribuicao menos tenue. (c) Distribuicao densa.

buicoes aleatorias de partıculas sao quase totalmente incoerentes, alem disso, em muitoscasos, tambem e possıvel considerar esse campo incoerente como uma funcao aleatoriaestacionaria no tempo, em determinados tamanhos de janela [16].

O princıpio para analise por ultra-som e a propagacao de pulsos no meio, cau-sando reflexoes devido as descontinuidades no caminho. Esses pulsos tem velocidadesrelativamente baixas, em torno de 1500 m/s para tecidos moles (Tabela 1.1), quandocomparadas com as ondas eletromagneticas, caso do raio-X, com valor igual a 3× 108

m/s. Devido a essa diferenca, os instrumentos eletronicos podem distinguir reflexoesde diferentes profundidades do corpo, pois o tempo consumido para o eco voltar ate oequipamento e da ordem de microsegundos, o que nao e possıvel com raio-X. Assim,pode-se reconstruir imagens com muito mais detalhes e sem nenhum dano (aparente)as celulas e tecidos.

Uma variedade de tecnicas sao usadas para inspecao, sendo as mais importantes:modos A, B, C, M e Doppler. O modo A mostra a amplitude do eco em funcao dotempo, o qual passa a ser proporcional a profundidade no meio. O modo B mostrauma imagem bidimensional, diferenciando os tecidos pelo brilho. O modo M e usadopara detectar movimento, basicamente das valvulas do coracao. O modo C e umatecnica que inspeciona profundidades constantes. Por ultimo, a inspecao por efeitoDoppler usa a variacao da frequencia da onda de ultra-som devido ao fluxo sanguıneo

6


para obtencao de informacoes [11, 50].Dessas tecnicas, as mais usadas sao o modo B e Doppler, devido a facilidade de

implementacao [22]. As imagens formadas sao apresentadas usando uma escala de tonsde cinza, sendo as propriedades diferentes do meio responsaveis pelas variacoes nasimagens. Pode-se citar a atenuacao e velocidade da onda de ultra-som e impedanciaacustica como algumas dessas propriedades. O modo A, apesar de menos usado, e oque apresenta maior quantidade de informacoes, pois analisa diretamente o sinal deRF (radiofrequencia), mas requer analises mais complexas.

Para as inspecoes, enviam-se pulsos de ultra-som de curta duracao para o tecidoe detectam-se os ecos que retornam das estruturas com diferentes impedancias. Pode-se usar o transdutor em modo pulso-eco, ou seja, como emissor e receptor, ou emuma configuracao com dois elementos piezeletricos. Quanto menor a duracao do pulso,maior a resolucao axial do sistema, isso pois a onda de ultra-som consegue interagircom estruturas menores. Um trem de pulsos de ultra-som tıpico, utilizado nas tecnicasde pulso-eco, e apresentado na Figura 1.2, [7].

Figura 1.2: Ilustracao de um trem de pulsos de ultra-som tıpico usado nas tecnicas depulso-eco. Apresentam-se a duracao do pulso, τ , e o perıodo de repeticao do pulso,PRP [7].

A duracao do pulso τ e bem curta, tipicamente menor que 1 µs, mas o perıodode repeticao do pulso PRP e relativamente longo, geralmente maior que 1 ms. Aintensidade e pressao instantaneas durante o pulso podem ser altas, algumas vezesalcancando dezenas de W/cm2 e alguns MPa, mas a intensidade media durante umciclo inteiro e baixa, geralmente menor que 100 mW/cm2 [7].

Essa duracao do pulso esta relacionada a largura de banda do transdutor. Peque-nos intervalos de tempo correspondem a grandes intervalos no espectro de frequencias,assim transdutores com largas bandas de frequencia sao preferıveis. Na fabricacao des-ses transdutores, conseguem-se bandas maiores para frequencias de operacao maiores,assim, essas duas caracterısticas estao geralmente atreladas.

Dessa forma, como o coeficiente de atenuacao aumenta linearmente com a frequenciapara tecidos moles, faz-se necessario usar um transdutor com frequencia de operacaomais baixa para se analisar uma estrutura mais profunda, em detrimento da resolucaoaxial. Para diminuir esse problema, alguns equipamentos possuem varios transdutores,escolhendo-se o mais apropriado para cada situacao, levando-se em conta a relacaoalcance-resolucao.

7


A resolucao lateral de um transdutor e proporcional ao produto da distancia focalpelo comprimento de onda na frequencia utilizada, dividido pela abertura (aperture).Portanto, aumentando-se a frequencia (diminuindo o comprimento de onda) ou aumen-tando a abertura, tem-se uma melhora na resolucao lateral. A Figura 1.3 (a) ilustra aresolucao axial e a Figura 1.3 (b), a lateral, considerando-se as imagens nomeadas como1 para um sistema com resolucoes suficientes e as 2 para um sistema com resolucoesinsuficientes [38].

Figura 1.3: Ilustracao das resolucoes axial e lateral de dois sistemas. As imagensnomeadas 1 sao de um sistema com resolucoes suficientes e as 2 de um sistema comresolucoes insuficientes. (a) Resolucao axial. (b) Resolucao lateral.

Um detalhe importante em imagens de ultra-som e a difracao nos tecidos, aqual causa divergencia da frente de onda. Pode-se somente focalizar o sistema emdeterminadas profundidades. Assim, devido a esse efeito, nao se consegue focalizar emregioes muito proximas ou muito distantes do transdutor [7].

Como exemplo do processo de obtencao de ecos de ultra-som de tecidos, apresenta-se a Figura 1.4(a) [38]. Neste, utiliza-se um transdutor que envia, pela aplicacao deum pulso de tensao eletrica em seus terminais, uma onda de ultra-som no tecido, aqual se propaga atraves de suas interfaces. Os ecos gerados por cada interface saodiferentes, devido as propriedades destas, como impedancia acustica, rugosidade e in-clinacao. Desconsiderando o espalhamento do tecido, o sinal de eco recebido em funcaoda distancia das interfaces e, de forma geral, como o apresentado na Figura 1.4(b).

Quando se considera o espalhamento, acrescentam-se informacoes sobre a consti-tuicao dos tecidos, alem de suas interfaces. Assim, os sinais de eco sao mais complexos,devido a interacao com partıculas do meio de tamanho comparavel ao comprimento deonda do ultra-som. A Figura 1.5 apresenta um exemplo de sinal de eco proveniente dainspecao do tecido da Figura 1.4, considerando-se tambem o espalhamento. As interfa-ces sao identificadas pelos picos e os valores intermediarios sao referentes a constituicaodo tecido.

Como o sinal de ultra-som e atenuado conforme avanca pelo meio, tecidos deigual poder de reflexao, mas a profundidades diferentes, gerarao ecos com intensidadesdiferentes, o que nao e desejado. Para se compensarem esses problemas, usa-se uma

8


(a)

(b)

Figura 1.4: Exemplo do processo de obtencao de ecos de ultra-som de tecidos. (a)Interacao da onda de ultra-som com as interfaces do tecido mostrando as diferencasnas reflexoes. (b) Sinal de eco recebido em funcao da distancia das interfaces, descon-siderando o espalhamento do tecido.

pre-amplificacao do sinal recebido com ganho variavel em funcao do tempo. Essa funcaoe geralmente logarıtmica, sendo mınima quando o pulso e emitido e progressivamentemaior com o avanco do tempo.

Inspecao Modo A

A inspecao modo A apresenta a amplitude do eco em funcao da distancia per-corrida pela onda de ultra-som, a qual e proporcional ao tempo, pois a velocidade doultra-som no meio e considerada constante. Essa tecnica foi desenvolvida em 1945 paradetectar falhas em metal, sendo, em 1950, utilizada pela primeira vez para examinartumores em tecidos humanos mortos, descobrindo-se variacoes de amplitudes dos ecosnos diversos tecidos [11].

O princıpio de funcionamento, como ja explicado, consiste em observar os ecosprovenientes das estruturas dos tecidos inspecionados, conforme exemplificado pelasFiguras 1.4 e 1.5. Como o eco de ultra-som consiste em um sinal com alguns ciclos,de amplitudes positivas e negativas, a informacao esta contida somente no envelopedo sinal, o qual e usado para obtencao das propriedades dos tecidos [11]. O tempo

9


Figura 1.5: Exemplo da Figura 1.4 considerando-se tambem o espalhamento do tecido.As interfaces sao identificadas pelos picos e os valores intermediarios sao referentes aconstituicao do meio.

que um eco demora para ser capturado pelo sistema, representa duas vezes a distanciado transdutor ate a estrutura do tecido que gerou esse eco. Portanto, sabendo-se avelocidade do ultra-som nesse meio, pode-se encontrar essa distancia.

A Figura 1.6(a) mostra um diagrama de blocos simplificado de um equipamentomodo A, conforme [9]. O elemento inicial do diagrama e o Gerador PRF (Pulse Re-petition Frequency), que cria um trem de pulsos com frequencia fixa, como mostra osinal A da Figura 1.6(b), o qual e usado como referencia para os outros componentes.

O Gerador de pulsos de tensao tem como funcao limitar a largura de pulso a seraplicada ao transdutor, sinal B da Figura 1.6(b). O componente Base de tempo cria umsinal dente-de-serra para varredura horizontal do mostrador CRT (Cathode Ray Tube),sinal C da Figura. O Pre-processamento do sinal realiza um condicionamento do sinalrecebido do transdutor, como mostra o sinal D, este constitui-se por um limitador detensao, um pre-amplificador e um gerador TGC (Time Gain Control), de controle deganho no tempo do pre-amplificador.

O bloco Demodulador e Amplificador e responsavel por fazer a deteccao da en-voltoria do sinal recebido e a sua suavizacao, sinal E, alem de fazer a compressaologarıtmica da faixa dinamica do mesmo para que este possa ser melhor observado peloser humano. O CRT recebe o sinal de varredura horizontal da base de tempo e devarredura vertical do Demodulador e Amplificador.

10


(a) (b)

Figura 1.6: Diagrama simplificado de um equipamento modo A. (a) Diagrama deblocos. (b) Exemplos dos principais sinais.

Inspecao Modo B

A inspecao modo B consiste em montar imagens que correspondam aos tecidosinspecionados, sendo a amplitude do eco representada pela intensidade ou nıvel decinza da imagem. Essa amplitude e, normalmente, comprimida logaritmicamente paraacentuar nıveis mais baixos. Ao inves de utilizar o sinal de eco para produzir umadeflexao vertical no monitor CRT, como no modo A, ele e usado para aumentar oudiminuir o brilho do feixe de eletrons, como mostra a Figura 1.7, identificando o objetopela diferenca no brilho da imagem. Com esse tipo de apresentacao, uma varreduralinear e aplicada ao eixo X do CRT, sendo o brilho do feixe ajustado para que se tornevisıvel somente quando um eco estiver presente. Dessa forma, a posicao dos pontosbrilhantes representam as distancias dos objetos [50].

Em alguns equipamentos Modo B atuais, a posicao do feixe de ultra-som tambeme monitorada, sendo apresentada junto com a informacao da amplitude dos ecos no mo-nitor CRT [7]. Esses equipamentos conseguem adquirir imagens rapido o suficiente paramonitorar o movimento dos orgaos. Alguns sistemas de imagens modernos possuemvarredura eletronica e focalizacao dinamica, o qual varre toda a regiao de interessecom o feixe focalizado, selecionando os melhores segmentos dos ecos para a geracao daimagem [23].

11


Figura 1.7: Exemplo de formacao de imagem modo B, utilizando-se os ecos para mo-dificar o brilho dos feixes de eletrons do monitor CRT.

A Figura 1.8 mostra um diagrama de blocos simplificado de um equipamentomodo B [7]. O elemento Condicionador do feixe e composto por uma matriz de atra-sos de fase e amplificadores, caso seja analogico, ou, se for digital, de arranjos depre-amplificadores e conversores A/D (analogico para digital). O TGC (Time GainControl) tem a mesma funcao que no aparelho modo A, ou seja, dar ganho conforme otempo, amplificando sinais provenientes de pontos mais distantes. O bloco de Proces-samento do sinal tem por funcao fazer a compressao da faixa dinamica, demodular efiltrar o sinal. O Conversor para CRT e o dispositivo de memoria digital que armazenaos dados que sao convertidos de seu formato padrao para um possıvel de ser exibidono CRT. Pode-se, antes de exibir o dado, fazer uma filtragem, como passa-banda oumapeamento em tons de cinza.

As imagens de ultra-som em modo B apresentam uma aparencia granular, cha-mada padrao de manchas (speckle pattern), causada por interferencias construtivas edestrutivas de pequenos ecos que chegam ao transdutor originadas do espalhamento emcomponentes do tecido. Como os pulsos de ultra-som, usados pelos aparelhos, sao ondasparcialmente coerentes [7], o padrao de manchas exibido pelo tecido possui informacoesrelevantes sobre sua microestrutura. A analise do speckle pode ter duas abordagensdiferentes: considera-lo como ruıdo ou como um sinal que contribui para a imagem[22]. No primeiro tipo, utilizam-se metodos para o eliminar, melhorando a qualidadeda imagem, e no segundo, analisa-se a sua distribuicao para buscar caracterısticas dotecido, sendo esta uma area promissora e pouco explorada [23].

Com base nos conceitos apresentados neste capıtulo introdutorio, este trabalhoaborda alguns metodos de processamento digital aplicados a sinais de ultra-som, inici-almente em sinais artificiais, criados a partir de modelos matematicos, e em seguida, asinais reais obtidos de phantoms que buscam simular propriedades de tecidos biologicos.O objetivo e obter alguns parametros que possam caracterizar o meio, pricipalmentequanto a sua regularidade.

Assim, no proximo capıtulo, uma revisao bibliografica apresentara alguns dosmetodos e modelos de tecidos biologicos utilizados normalmente, mostrando tambemos parametros dos meios que caracterizam sua regularidade e qual a importancia dessapropriedade para a inspecao por ultra-som. Serao apresentados tambem alguns exem-plos simples da aplicacao dos metodos para uma melhor definicao do problema.

No terceiro capıtulo, explica-se em detalhes os conceitos teoricos de processa-mento de sinais e filtragem adaptativa empregados nas analises dos sinais. Inicia-se

12


Figura 1.8: Diagrama de blocos simplificado de um equipamento modo B.

com uma introducao a processos estocasticos, que sao a base para o resto dos estudos.Em seguida, apresentam-se os Modelos Autorregressivos (AR), os quais foram usadosna parametrizacao dos sinais de ultra-som para extracao das caracterısticas de regula-ridade do meio. Alguns conceitos de filtragem adaptativa, como filtros de Wiener, depredicao linear e LMS tambem sao apresentados, buscando uma ampla descricao dateoria envolvida. Apos isso, explica-se o metodo de rastreamento utilizado nos polos dadecomposicao AR para observacao das mudancas quando se variavam os parametro domeio. Por ultimo, apresenta-se o metodo de simulacao de Monte Carlo utilizado parase obter resultados estatısticos das hipoteses levantadas e os aspectos computacionaisnecessarios para o desenvolvimento dos algoritmos.

No quarto capıtulo, apresentam-se os materiais e metodos utilizados no desen-volvimento deste trabalho. Primeiramente, explicam-se como foram feitos os meios deteste utilizados para obtencao de sinais reais. Em seguida, os equipamentos e softwa-res de aquisicao de dados. Apos estes, apresentam-se os softwares de analise, com osmetodos numericos e de programacao utilizados, junto com as ferramentas emprega-das. Por ultimo, apresentam-se as simulacoes, com a criacao dos ecos artificiais pelosmodelos matematicos e os metodos de analise.

No quinto capıtulo, apresentam-se os resultados obtidos neste estudo. Primeira-mente, aqueles obtidos na simulacao dos sinais artificiais; em seguida, os resultados dasimulacao de Monte Carlo, mostrando a sensibilidade do sistema a variacoes em dois

13


parametros dos meios: a regularidade e a relacao sinal-ruıdo. Por ultimo, apresentam-se os resultados das analises aplicadas a sinais de ultra-som reais obtidos de phantoms.No sexto capıtulo, apresentam-se as conclusoes e discussoes dos resultados apresen-tados no capıtulo anterior, juntamente com as sugestoes de proximos passos para acontinuacao do trabalho.

14

Capıtulo 2

Revisao Bibliografica

Neste capıtulo, serao apresentadas revisoes bibliograficas sobre o modelamento de te-cidos biologicos, a caracterizacao de tecidos por espacamento medio dos espalhadores,MSS (Mean Scatterer Spacing), e os metodos para se encontrar esse parametro.

O modelamento apresentado na Secao 2.1 considera tecidos ditos como moles, quee o caso de fıgado e baco por exemplo. Esses tecidos podem ser modelados por umasuspensao de partıculas, apresentando estruturas regulares com propriedades bem defi-nidas [22]. Segundo [18], os tecidos de fıgado podem apresentar dois tipos de doencas:localizadas ou difusas. A primeira se caracteriza pela concentracao em pequenos pon-tos em um ou nos dois lobulos do orgao, sendo o resto do tecido normal. Na segunda,pelo menos um dos lobulos sao afetados por inteiro, sendo a doenca espalhada. Essaclassificacao nao faz distincao por importancia ou grau de avanco da doenca.

Assim, como esses tipos de tecido apresentam estruturas regulares distribuıdashomogeneamente, variacoes causadas por doencas difusas podem ser detectadas pormudancas em parametros quantitativos obtidos de sinais retroespalhados de ultra-som.Esse e o caso do MSS, que tende a apresentar um aumento em tecidos nessas condicoes[4, 21, 22, 29, 35], como explicado na Secao 2.2.

2.1 Modelamento de Tecidos Biologicos

A apresentacao de um modelo de tecidos biologicos que ajude na descricao da formacaodos ecos de ultra-som retroespalhados e de grande importancia para o estudo dessessinais e a determinacao de parametros quantitativos para a caracterizacao dos tecidos.Isso pois, geralmente, dados coletados por ultra-som nao podem ser interpretados demaneira quantitativa sem estarem embasados em um modelo. Segundo [1], o objetivoprincipal de um modelo e prever a resposta de um sistema as propriedades do meio eanomalias em sua estrutura. Tambem e necessario incluir a configuracao do sistema deaquisicao de sinais, considerando a resposta do transdutor e as interacoes da onda deultra-som.

A interacao entre o sinal de ultra-som e as estruturas espalhadoras dos tecidosbiologicos sao analisadas de duas formas diferentes: primeira, considerando-se um mo-delo determinıstico para descrever as estruturas regulares, ou segunda, com um modeloestatıstico, o qual descreve o tecido como um conjunto de espalhadores distribuıdos deforma aleatoria. Entretanto, nenhum deles modela completamente os tecidos, que, de

15

Capıtulo 2 - Revisao Bibliografica

forma mais realista, deveriam ser distribuicoes aleatorias com certo grau de regulari-dade estatıstica [21].

Em [21], um modelo e proposto objetivando uma melhor compreensao dos pro-cessos responsaveis pela formacao dos ecos. Considera-se que, apesar dos tecidos seremtridimensionais, de forma simplificada, eles podem ser descritos por um modelo unidi-mensional. Entao, fazendo-se a inspecao do meio com um feixe estreito de ultra-som,pode-se dizer que o sinal e espalhado por uma linha de espalhadores, com posicao eforca aleatorias.

Dessa forma, para se obter o sinal retroespalhado, considera-se o tecido como umacolecao de espalhadores aleatoriamente distribuıdos, filtrados por sistema invarianteno tempo que representa o equipamento de medida e as propriedades de espalhamentomedias do meio. Assim, a resposta impulsiva do tecido, h(t), e definida pela Equacao2.1, considerando-se um sistema linear, o que e valido, pelo menos, para o caso deanalise de pequenos sinais [21].

h(t) =∑

i

gi(t− ti) (2.1)

Nessa equacao, gi(t) e uma forma de onda selecionada de um grupo de formas deonda fisiologicamente restrito, que leva em conta a geometria e a forca de espalhamentodo i-esimo espalhador, sendo ti, o atraso de tempo da resposta, o qual depende daposicao aleatoria da partıcula. Considerando que a resposta impulsiva seja um processoWSS, e que as formas de onda gi(t) sejam independentes da sequencia de tempo ti,chega-se a Equacao 2.2.

Egh(t) = g(t) ∗∑

i

δ(t− ti) (2.2)

O sımbolo ∗ representa o operador convolucao e Eg. significa a esperanca doconjunto de ondas gi(t). O somatorio dos impulsos deslocados no tempo pode serinterpretado como uma funcao amostrada de um processo de entrada, reproduzindouma configuracao particular da distribuicao espacial dos espalhadores. Assim, a res-posta do sistema e facilmente obtida pela convolucao desse somatorio com g(t), comoilustrado pela Figura 2.1.

A modelagem usando funcoes delta e uma simplificacao matematica util quando ointeresse esta somente no tempo de chegada do eco. Entao, como o objetivo da analisee a compreensao dos efeitos da regularidade nas distancias entre os espalhadores, aspropriedades estatısticas da variavel aleatoria (conceito a ser definido no Capıtulo 3)τ = ti − ti−1, chamada de tempo de voo, a qual representa o intervalo entre a chegadade dois ecos consecutivos, devem ser analisadas.

A regularidade do modelo quer dizer que o desvio padrao da variavel τ e pequenoquando comparado com o valor medio, τ , sendo o modelo mais regular, quanto menorfor essa relacao. Essa regularidade e caracterizada por um pico estreito na funcao dedensidade de probabilidade (Capıtulo 3) de τ ao redor de seu valor medio [21].

Para o caso de baixa regularidade, espera-se que essa funcao se aproxime deuma distribuicao de Poisson com forte caracterıstica aleatoria. Para alta regularidade,a funcao de densidade de probabilidade pode ser aproximada por uma distribuicao

16

2.1 Modelamento de Tecidos Biologicos

Figura 2.1: Ilustracao do modelo unidimensional apresentado em [21] considerando ostecidos biologicos como uma distribuicao aleatoria de espalhadores (scatterers) comforma e forca de espalhamento aleatorios.

Gamma, que especifica intervalos τ independentes e distribuıdos igualmente, comomostra a Equacao 2.3 [31].

f(τ) =1

Γ(n)(τ /n)nτn−1e−nτ/τ , τ ≥ 0 (2.3)

Γ(n) =

∫ ∞

0

tn−1e−tdt = (n− 1)!

Como se pode ver, a distribuicao Gamma depende de duas constantes: o valormedio do intervalo de chegada, τ , e um ındice n, que modifica a variancia de τ , poisσ2

τ = τ 2/n. Essa variancia mede o grau de irregularidade do processo, para n = 1,o processo e bem proximo de um Poisson aleatorio e, para valores grandes de n, oprocesso e bem regular.

Em [17], considera-se tambem um modelo de tecido composto por uma distri-buicao de espalhadores aleatorios com aplicacoes em tecidos de testıculo e prostata.Nesse trabalho, apos correcoes de influencias dependentes do sistema, utilizam-se analisesespectrais dos sinais de RF usados na formacao de imagens modo B para acrescentarinformacoes sobre os tecidos, alem das obtidas pelas imagens.

Pelo uso do modelo, realizaram-se testes em phantoms de agar e po de grafitepara se obter as funcoes de correcao. Os parametros atenuacao e coeficiente de retro-espalhamento sao estimados como valores medios em regioes de interesse dos tecidos.Fazem-se as analises espectrais dentro da largura de banda do transdutor usado, pois,esses parametros sao dependentes da frequencia. As estimativas sao satisfatorias, re-sultando em um aumento das probabilidades de deteccao de doencas nos tecidos.

Nos trabalhos [4, 29, 35], uma abordagem bem parecida e utilizada, mas que apre-senta algumas vantagens. O estudo e feito para tecidos de fıgado, buscando diferenciaros saudaveis dos cirroticos. Para isso, o modelo usado considera dois tipos de estru-turas espalhadoras: a primeira, de espalhadores difusos, os quais possuem distribuicaoaleatoria, caracterizando a parte difusa do tecido; a segunda, espalhadores regulares,distribuıdos de maneira periodica ou quase-periodica. Da mesma forma que em [21],o eco y(t) e modelado como a convolucao da resposta do tecido, h(t) com o pulso de

17


ultra-som, p(t), como apresentado na Equacao 2.4, com uma diferenca na distribuicaodos espalhadores, Equacao 2.5.

y(t) = h(t) ∗ p(t) (2.4)

h(t) =Nr∑n=1

rn(t− trn) +

Nd∑n=1

dn(t− tdn) (2.5)

Nessa equacao, Nd e o numero de espalhadores difusos, Nr, o numero de espa-lhadores periodicos, trn e tdn sao os atrasos de tempo dos espalhadores ate o receptore , por ultimo, rn(t − trn) e dn(t − tdn) sao as forcas relativas dos espalhadores. Aposicao dos espalhadores difusos e dada por uma distribuicao uniforme e dos periodicose Gamma. As forcas de espalhamento sao distribuicoes uniformes. Em [4], as analisesdos ecos foram feitas usando decomposicoes WOLD, que apresentaram bons resultados,sendo possıvel diferenciar os fıgados saudaveis de outros com doencas pelo parametrode espacamento medio dos espalhadores, o qual sera apresentado na proxima secao.

(a) (b)

Figura 2.2: (a) Esquema da constituicao do tecido de fıgado humano. (b) Corte de umfıgado humano mostrando as estruturas que o compoem (Figuras retiradas de [38]).

Segundo [29], algumas consideracoes precisam ser feitas para que o modelo apre-sente bons resultados: primeiro, as forcas relativas do espalhadores, rn e dn, para cadaespalhador, devem ser variaveis aleatorias (VA) independentes; segundo, os atrasosde tempo dos espalhadores difusos, tdn, tambem devem ser VA independentes, dis-tribuıdas uniformemente em 0, 2π; terceiro, os espalhadores periodicos devem serdescorrelacionados com os difusos e, por ultimo, o eco, y(t), precisa ser um processoWSS.

Como ja dito, o fıgado e considerado um tecido mole, o qual e organizado empequenas unidades histologicas chamadas de lobulos hepaticos, que podem ser relaci-onados aos componentes periodicos do modelo proposto [35]. Esses lobulos possuem

18

2.2 Espacamento Medio dos Espalhadores

uma parede de tecido conjuntivo, em formato pentagonal ou hexagonal, com arestasdadas pelos espacos porta, como mostra o esquema da Figura 2.2(a) (retirada de [38]),sendo estas estruturas espalhadas periodicamente pelo tecido [38]. As outras estruturaspresentes sao associadas aos espalhadores difusos. A Figura 2.2(b) (retirada de [38])apresenta um corte de um fıgado humano mostrando essas estruturas.


Conforme apresentado anteriormente, tecidos moles podem ser modelados, do ponto devista do eco retroespalhado de ultra-som, como um conjunto de dois tipos de espalha-dores, um distribuıdo (quase) periodicamente e outro, aleatoriamente. Para o caso dostecidos de fıgado, a parte regular e devida a estruturas chamadas lobulos hepaticos, asquais sao distribuıdas de forma homogenea no meio.

Figura 2.3: Ilustracao do espacamento medio dos espalhadores, considerando o modelounidimensional apresentado anteriormente.

Assim, na inspecao desses tecidos por ultra-som, pode-se retirar, do sinal de eco,um parametro quantitativo do meio chamado espacamento medio dos espalhadores,MSS (Mean Scatterer Spacing), o qual relaciona o tempo medio de chegada ao receptor(transdutor) dos ecos periodicos com o espacamento das estruturas regulares, comoilustrado pela Figura 2.3. Com base no modelo utilizado, tem-se que o MSS e definidopelos ecos dos espalhadores periodicos.

Segundo o modelo proposto por [21], o valor medio do intervalo de chegada dedois ecos consecutivos dos espalhadores regulares e dado por τ , assim, considerando-sea velocidade do ultra-som no meio como c, tem-se o valor do MSS, conforme a Equacao2.6.

MSS =c τ

2(2.6)

Em varios trabalhos [4, 21, 29, 35, 44], bons resultados foram obtidos usando-se oMSS para caracterizar tecidos e doencas difusas. Nesses estudos, a forma mais comum

19


para estimativa desse parametro e pelo uso do espectro do eco, usando a amplitude oua fase. Mas, como o eco e a convolucao do pulso com o modelo do meio, o espectrose apresenta modificado pela resposta em frequencia do transdutor, alem de possuirinterferencias da parte difusa do tecido.

A Figura 2.4 mostra dois casos de resposta impulsiva do tecido e seus espectrosde potencia, segundo o modelo apresentado. A diferenca entre eles esta no desviopadrao do tempo de voo, τ , que mede a regularidade do meio. Quando este parametro,tambem chamado de jitter [22], e pequeno em relacao a media, tem-se que o meio equase-periodico, mas conforme ele aumenta, perde-se essa caracterıstica. O valor dojitter e dado como uma porcentagem da media, assim, na figura, usam-se 5 e e 50%.

Como se pode ver nessas figuras, a caracterıstica de quase-periodicidade do meioimplica em raias no espectro. Isso, pois a parte regular pode ser considerada como umtrem de impulsos, o qual possui um espectro tambem como um trem de impulsos [30].Entao, se o modelo fosse totalmente periodico, o espectro tambem o seria. Mas, devidoao jitter, os picos de alta frequencia do espectro sao bem atenuados, e, quando essevalor e alto, o espectro possui praticamente so o primeiro pico, sendo os harmonicosconfundidos com os outros elementos do sinal.

Dessa forma, se fosse possıvel utilizar o espectro da resposta do meio para seencontrar o valor do MSS, seria necessario achar o valor de frequencia correspondenteao primeiro pico do espectro, sendo o espacamento medio dado pela Equacao 2.7 [22].

MSS =c

2 fMSS(2.7)

Entretanto, na formacao do eco, leva-se em consideracao a resposta impulsiva dotransdutor, assim, o espectro de potencia do eco, Sy(w) tambem possui componentesda resposta em frequencia do pulso, P (w), como mostra a Equacao 2.8, sendo Sh(w), oespectro do modelo do meio [29]. A Figura 2.5 (a) apresenta o espectro de potencia daresposta do transdutor e (b), o espectro do eco, mostrando o deslocamento provocadopela frequencia de ressonancia e largura de banda do transdutor.

y(t) = h(t) ∗ p(t) ⇒ Sy(w) = |P (w)|2Sh(w) (2.8)

Para a estimativa do MSS a partir do eco, varios metodos podem ser usados.Alguns trabalhos obtiveram bons resultados usando metodos diferentes, os quais podemser classificados em dois grupos [35]: os que usam a informacao contida na amplitudedo espectro do eco e os que usam tanto a amplitude como a fase.

No primeiro grupo, pode-se citar Landini e Verrazzani [21] que usaram o ceps-trum da magnitude do espectro para converter o efeito multiplicativo do espectro daresposta do transdutor em um efeito aditivo. Wear et al., [48, 49] usaram o modeloautorregressivo pelo metodo de Burg para estimar o espectro do eco, obtendo bonsresultados para colecoes de dados de pequenos comprimentos usando phantoms feitosde agar e esferas de vidro e com amostras de fıgado humano.

20


(a)

(b)

Figura 2.4: Exemplos da resposta impulsiva do tecido usando o modelo, junto comseus espectros de potencia, para jitter de (a) 5% e (b) 50%. A variacao da amplitudedos ecos, nestes exemplos, e imperceptıvel visualmente.

21


No outro grupo, Varghese e Donohue [46] utilizaram a tecnica de autocorrelacaoespectral, que tambem apresenta bons resultados para colecoes de dados de pequenocomprimento, aperfeicoando essa funcao com suavizamento em frequencia em um tra-balho seguinte [45]. A desvantagem desse metodo e o grande custo computacionaldevido a transformacoes bidimensionais para frequencia, o que inviabiliza as aplicacoesem tempo real.

Figura 2.5: (a) Espectro de potencia da resposta do transdutor. (b) Espectro depotencia do eco mostrando o deslocamento provocado pela frequencia de ressonancia elargura de banda do transdutor.

Outro metodo que usa a informacao de fase do espectro e apresentado por Simonet al. [35], o qual e mais simples que o anterior por usar transformacoes unidimensio-nais. Essa tecnica, que tambem e utilizada por Maciel [22], consiste em realizar umatransformacao quadratica do eco antes de se encontrar o espectro, de forma a acentuara primeira harmonica, tornando-a mais facil de ser detectada.

Por sua simplicidade de implementacao, esse metodo foi adotado para ser uti-lizado neste estudo. O espectro de potencia do sinal de eco quadrado, y2(t), in-dica a correlacao entre os componentes espectrais de mesma frequencia, entao, se aspartıculas estao espacadas irregularmente, os valores esperados de correlacao cruzadaem frequencia sao insignificantes [22]. Assim, o primeiro pico do espectro que repre-senta a frequencia fundamental dos trem de impulsos dos componentes regulares serao que apresentara maior correlacao, destacando-se no espectro, podendo ser detectadocomo maximo do sinal.

22


Figura 2.6: (a) Espectro de potencia do exemplo de sinal de eco apresentado na fi-gura anterior. (b) Espectro de potencia do sinal ao quadrado, mostrando o primeiroharmonico com maior amplitude, correspondendo a fMSS, e a diminuicao do efeito daresposta do transdutor.

A Figura 2.6 (a) apresenta o espectro de potencia do exemplo de eco apresentadona Figura 2.5 (b) para ser comparado com o espectro do eco ao quadrado, (b). Pode-sever que o primeiro harmonico possui a maior amplitude, correspondendo ‘a frequenciade MSS, fMSS, e que o efeito da resposta do transdutor foi diminuıdo.

23

24

Capıtulo 3

Teoria

Os sinais utilizados em estudos de sistemas podem ser classificados, de modo geral,em dois tipos: determinısticos ou aleatorios. Um sinal determinıstico e aquele quepode ser exatamente reproduzido em diversas medidas, podendo ser expresso por al-guma formula matematica explıcita. Em oposicao, para um sinal aleatorio isso nao epossıvel, este tipo so pode ser descrito probabilisticamente ou utilizando suas proprie-dades estatısticas, como seus momentos [13].

Poucos sinais utilizados na pratica podem ser considerados determinısticos, sendoexemplos a saıda de fontes de energia, osciladores ou geradores de sinais, os quais pre-cisam apresentar sinais sempre iguais, isso, desconsiderando-se os erros e variacoes nor-mais dos equipamentos. A maioria dos sinais usados sao classificados como aleatorios,ou como alguns autores denominam [41] sinais portadores de informacao, pois elesnao poderiam ser chamados verdadeiramente de aleatorios, sendo tratados por suascaracterısticas medias.

3.1 Processos Estocasticos

Quando se tratam de sinais aleatorios, um conceito importante e o de variavel aleatoria(VA). Para uma explicacao consistente, algumas definicoes precisam ser feitas. A pri-meira e a de espaco amostral, Ω, o qual se caracteriza por todos os resultados possıveisde um experimento, como, por exemplo, o lancamento de um dado, que possui espacoamostral de seis eventos elementares ω1, ω2, ω3, ω4, ω5, ω6, correspondendo a cada umdos possıveis lados que podem ficar para cima quando o dado se estabilizar. A segundadefinicao e a probabilidade de um evento, a qual e a chance de um determinado eventoocorrer em um experimento, no caso do dado, a probabilidade do lado 1 e denotadopor Pω1 = 1/6.

Entao, de posse desses conceitos, pode-se definir variavel aleatoria como o pro-cesso de se designar um numero x(ω)para cada evento ω, [31]. No exemplo anterior,uma possıvel VA seria criada designando-se os numeros de um a seis para os eventos ω1 aω6, assim, uma forma de se apresentar a probabilidade do lado 1 seria Px = 1 = 1/6.

O exemplo anterior apresenta uma VA discreta mas, uma visao mais geral baseia-se no tipo contınuo, sendo o espaco amostral Ω mapeado no eixo real (ou complexo).Um exemplo desse caso e um sinal de ruıdo eletrico, que possui amplitude real e

25

Capıtulo 3 - Teoria

aleatoria. Nas aplicacoes de processamento de sinais, uma descricao probabilısticade VA e mais interessante do que uma caracterizacao estatıstica de eventos [13].

Duas formas de se descrever uma VA sao as funcao de distribuicao e de densidadede probabilidade, F (x) e f(x) respectivamente. A primeira representa a probabilidadeda variavel x ter um valor menor que x, Px ≤ x, e a segunda e a derivada dessafuncao, dF (x)/dx. Porem, em muitos casos, e difıcil se trabalhar com a funcao dedensidade (ou distribuicao) de probabilidade por esta ser desconhecida. Nesses casos,e possıvel utilizar somente os momentos da distribuicao, principalmente os de primeirae segunda ordem [41].

O momentos mais usados sao o de primeira ordem, conhecido como media, e o mo-mento central de segunda ordem, conhecido como variancia, definidos pelas Equacoes3.1 e 3.2 respectivamente, sendo E. o operador de valor esperado.

µ = Ex =

∫ ∞

−∞xf(x)dx (3.1)

σ2 = E(x− µ)2 =

∫ ∞

−∞(x− µ)2f(x)dx (3.2)

Outras caracterısticas estatısticas usadas quando se lidam com duas VA, x e y porexemplo, sao os momentos conjuntos (joint moments). Os principais sao a correlacaoe a covariancia, definidos nas Equacoes 3.3 e 3.4 respectivamente, sendo µx e µy asmedias das VA e o operador .∗ equivalente a conjugado complexo.

rxy = Exy∗ (3.3)

cxy = E(x− µx)(y − µy)∗ (3.4)

E comum, utilizar-se um valor de covariancia normalizado, chamado coeficientede correlacao, definido pela Equacao 3.5.

ρxy =cxy

σxσy

(3.5)

Um tipo importante de VA e a gaussiana, que possui funcao de densidade deprobabilidade da forma da Equacao 3.6, sendo totalmente caracterizada pela media evariancia da variavel, µ e σ2 [31].

f(x) =1√

2πσ2exp

−(x− µ)2

2σ2

(3.6)

Um processo aleatorio ou estocastico pode ser descrito como a evolucao temporal(definida em um intervalo de observacao) de um fenomeno estatıstico de acordo comleis probabilısticas, portanto, sem o conhecimento exato de sua evolucao [14]. Algunsexemplos de processos aleatorios sao: sinais de ultra-som, voz, telecomunicacoes eruıdos.

Uma definicao mais teorica desses processos pode ser apresentada como um ma-peamento do espaco Ω em um conjunto de sinais discretos no tempo, ~x(n), sendoentao, uma sequencia indexada de variaveis aleatorias [13]. Assim, para cada eventoωi do espaco Ω, existe um sinal discreto correspondente no tempo, xi(n), chamado de

26


serie temporal, o qual e composto por varios resultados diferentes do mesmo eventoωi. Pode-se pensar que, para um dado n = k, o valor do processo, ~x(n) e uma variavelaleatoria, a qual pode ser chamada xk, que possui uma designacao para cada eventoωi, nominalmente xk(i), conforme ilustrado na Figura 3.1.

Figura 3.1: Ilustracao da formacao de um processo estocastico, ~x(n), considerando-seo espaco amostral Ω, cujos eventos ωi geram as variaveis aleatorias xk(i).

Conforme apresentado anteriormente, cada uma dessas VA possuem uma funcaodensidade (e distribuicao) de probabilidade. Entao, para se caracterizar completa-mente um processo, precisa-se das funcoes de densidade (ou distribuicao) conjuntapara qualquer colecao de VA [31].

Como um processo estocastico e uma colecao de VA, pode-se calcular a sequenciade valores medios e variancia de cada VA, as quais sao dependentes do tempo n, comomostram as Equacoes 3.7.

µ(n) = Ex(n)

σ2(n) = E|x(n)− µ(n)|2 (3.7)

Pode-se encontrar tambem a autocorrelacao e autocovariancia do processo, deforma similar ao apresentado anteriormente (Equacoes 3.3 e 3.4), considerando-se duasVA quaisquer do processo, x(n) e x(k), como mostram as Equacoes 3.8.

r(n, k) = Ex(n)x∗(k)

c(n, k) = E[x(n)− µ(n)][x(k)− µ(k)]∗ (3.8)

27


No caso de processos com media nula, as duas equacoes se tornam iguais, oque e normalmente considerado neste trabalho. Quando n = k, a autocorrelacao setorna igual a variancia. A correlacao cruzada de dois processos pode ser encontradasubstituindo-se na equacao anterior uma das VA por uma do outro processo.

A forma do processo estocastico utilizado neste trabalho caracteriza-se por dis-creto e igualmente espacado no tempo. Uma propriedade importante e a estacionarie-dade estrita, enunciada como a invariancia das caracterısticas estatısticas a deslocamen-tos temporais. Esse caso ideal e muito raro na pratica, sendo geralmente substituıdopelo conceito de Estacionariedade em Sentido Amplo (WSS - Wide Sense Stationarity).

Um processo e dito WSS quando as seguintes tres condicoes sao satisfeitas: pri-meira, a media do processo e constante; segunda, a autocorrelacao r(n, k) so dependeda diferenca n − k; e, terceira, a variancia do processo e finita [31]. Processos WSSpossuem tres propriedades interessantes [13]:

• A sequencia de autocorrelacao do processo e uma funcao conjugada simetrica dek, r(k) = r(−k)∗.

• A autocorrelacao do processo para um atraso k = 0 e igual ao valor medioquadratico do processo: r(0) = E|~x(n)|2.

• A magnitude da autocorrelacao de um processo WSS para um atraso k e limitadapelo valor para k = 0, r(0) ≥ |r(k)|.

Entretanto, os calculos dos momentos sao definidos usando funcoes de densidadede probabilidade (conforme as Equacoes 3.1 e 3.2). Mas, mesmo quando essas funcoesnao sao conhecidas, pode-se estimar esses momentos. Para isso, utiliza-se uma estima-tiva nao influenciada (unbiased estimator) [5, 13, 14].

O operador valor esperado nas formulacoes apresentadas de processos estocasticose considerado como media de grupo (ensemble average), realizado sobre as VA quecompoem um processo. Pode-se, tambem, realizar a media ao longo do tempo (timeaverage), a qual converge (em sentido medio quadratico para processos WSS [14]) paraa media de grupo, conforme o teorema ergodico da media. Assim, pode-se calcular osmomentos de uma forma alternativa, como apresentado para a media e autocorrelacaonas Equacoes 3.9 e 3.10, sendo N o numero total de amostras usadas na estimativa ex(n) uma serie temporal.

µ =1

N

N−1∑n=0

x(n) (3.9)

r(k) =1

N − 1

N−1∑n=0

x(n)x(n− k), 0 ≤ k ≤ N − 1 (3.10)

Um tipo de processo estocastico discreto no tempo muito importante em proces-samento de sinais e o ruıdo branco. Considera-se um processo WSS, ~v(n) como ruıdobranco, se a funcao de autocovariancia for zero para todos os k 6= 0, como mostra aEquacao 3.11. Assim, esse processo se caracteriza por ser uma simples sequencia deVA nao correlacionadas, as quais possuem variancia σ2

v . Como essa definicao so se ba-seia na forma do momento de segunda ordem, existem varios tipos de ruıdos brancos,

28


sendo um dos mais comuns, o ruıdo branco gaussiano, cujas VA possuem distribuicaogaussiana.

c(k) = σ2vδ(k) (3.11)

Uma forma comum de se apresentar as sequencias de autocorrelacao de um pro-cesso estocastico discreto no tempo e a matriz de correlacao. Considerando-se umprocesso WSS composto pela serie temporal x(n), apresentada na Equacao 3.12, esuas versoes atrasadas, obtem-se a matriz de correlacao, conforme a Equacao 3.13,sendo .T e .H os operadores de transposicoes matricial e hermitiana (transpostoconjugado) respectivamente.

~x(n) = [x(n), x(n− 1), . . . , x(n− p− 1)]T (3.12)

~R = E~x(n)~xH(n) =

r(0) r(1) . . . r(p− 1)

r(−1) r(0) . . . r(p− 2)...

.... . .

...r(−p + 1) r(−p + 2) . . . r(0)

(3.13)

Devido a importancia dessa matriz, algumas de suas propriedades sao apresenta-das para o caso de um processo estocastico WSS discreto no tempo [14]. A matriz decorrelacao e:

• Hermitiana, ou seja, e igual a sua transposta conjugada.

• Toeplitz, ou seja, os elementos ao longo de cada diagonal (a principal e as para-lelas) sao iguais.

• Nao-negativa definida e quase sempre positiva definida.

• Nao-singular, ou seja, seu determinante nao e zero, isso devido a presenca deruıdo aditivo.

A analise de Fourier e uma importante ferramenta tanto para sinais determinısticoscomo para aleatorios, normalmente encontrada em ciencias fısicas e engenharias. Mas,como esses sinais aleatorios sao um conjunto de sinais discretos, a transformada deFourier deve ser expressa em termos de uma media de grupo (ensemble) [13]. Conside-rando que a sequencia de autocorrelacao de um processo WSS e determinıstica, pode-secalcular sua transformada. As Equacoes 3.14 e 3.15 apresentam, respectivamente, oespectro de potencia de um processo ~x(n) (junto com a transformada z) e sua inversa.

Px(ejw) =

∞∑k=−∞

r(k)e−jwk (3.14)

Px(z) =∞∑

k=−∞

r(k)z−k

29


r(k) =1

2π

∫ π

−π

Px(ejw)ejwkdw (3.15)

Algumas propriedades interessantes do espectro de potencia de um processo WSSsao apresentadas a seguir [13]:

• O espectro de potencia e uma funcao real de w, pois a autorrelacao e simetricaconjugada.

• O espectro sera par, Px(ejw) = Px(e

−jw), se o processo ~x(n) for real.

• O espectro de potencia e nao-negativo.

• A potencia de um processo de media nula e proporcional a area sob a curva doespetro de potencia:

E|~x(n)|2 =1

2π

∫ π

−π

Px(ejw)dw

• Os autovalores da matriz de correlacao de um processo com media nula, λi, saolimitados pelo valores maximo e mınimo do espectro de potencia:

minwPx(ejw) ≤ λi ≤ maxwPx(e

jw)

3.2 Modelos Autorregressivos

A estatıstica de segunda ordem de um processo estocastico pode ser apresentada al-ternativamente por suas funcoes de autocorrelacao ou densidade espectral de potencia(PSD), as quais sao descricoes nao-parametricas [19]. Mas, formas parametricas podemser usadas para caracterizar esses processos, como e o caso dos modelos autorregressivos(AR).

Pelo uso de modelos parametricos em processos aleatorios pode-se, geralmente,obter uma melhor estimativa espectral do processo, estimando-se os parametros domodelo a partir das amostras disponıveis de dados. Estes metodos apresentam umaresolucao aparentemente maior para estimativas espectrais do que tecnicas classicas,quando se escolhe um modelo apropriado ao processo [19]. A vantagem e poder fazersuposicoes mais realistas a respeito do processo fora do intervalo de observacao sem fa-zer uso de funcoes de janelas. Alem disso, pode-se aplicar esses metodos para processosvariantes no tempo, se essa variacao for lenta [37].

O termo modelo e usado para qualquer hipotese que possa ser aplicada para des-crever as leis ocultas que, presumidamente, governam a geracao de dados fısicos deinteresse [14]. Modelos parametricos muito usados para aproximar processos aleatoriossao as funcoes de transferencia racional, nas quais uma sequencia de saıda y[n], alta-mente correlacionada e que modela as amostras de dados, e relacionada com uma en-trada x[n], normalmente ruıdo branco, por uma equacao linear de diferencas, Equacao3.16.

30


y[n] = −p∑

k=1

a[k] y[n− k] +

q∑k=0

b[k] x[n− k] (3.16)

Esse tipo de modelo e conhecido como autorregressivo media movel, ARMA (Au-toregressive Moving Average), conforme apresentado na Figura 3.2. Esse modelo e cha-mado racional devido ao formato da funcao do sistema H(z), como mostra a Equacao3.17.

H(z) =B(z)

A(z)=

q∑k=0

b[k] z−k

p∑k=0

a[k] z−k

(3.17)

Uma consideracao importante e que todos os polos dessa funcao estejam dentrodo cırculo unitario do plano z, para que H(z) seja estavel e causal [30]. Esses modelostambem possuem a propriedade de serem lineares e invariantes no tempo, pois assim,y[n] e uma descricao valida para um processo WSS [19].

Figura 3.2: Diagrama de um modelo ARMA para um processo estocastico, consi-derando a entrada x[n], tipicamente ruıdo branco, e a saıda y[n]. Para ilustracao,considera-se p = q.

A entrada, x[n], como ja apresentado, e um processo aleatorio de ruıdo brancogeralmente gaussiano com media nula e funcao de autocorrelacao dada pela Equacao3.18. Normalmente, os coeficientes a[0] e b[0] sao considerados iguais a 1, porquequalquer ganho do filtro pode ser incorporado no valor de σ2

x.

rx(n, k) = E[x(n) x∗(k)] =

σ2

x, k = n0, caso contrario

(3.18)

31


Considerando a Equacao 3.17, pode-se encontrar o espectro de potencia para asequencia de saıda do filtro ARMA, Py(z), como apresentado na Equacao 3.19. Comoa entrada e um ruıdo branco, o seu espectro e Px(z) = σ2

x. Avaliando-se o espectrono cırculo unitario, z = ejw, encontra-se a funcao de densidade espectral de potencia,Py(e

jw), apresentada na Equacao 3.20. Assim, o processo y[n] e chamado ARMA deordem (p, q), possuındo um PSD com 2p polos e 2q zeros.

Py(z) = H(z) H∗(1/z∗)Px(z) (3.19)

Py(ejw) = σ2

x

|B(ejw)|2

|A(ejw)|2(3.20)

Dois casos especiais de processos ARMA sao: media movel, MA (Moving Ave-rage), e autorregressivo, AR (Autoregressive). O primeiro surge quando os coeficientesa[k] sao iguais a zero, exceto a[0] = 1, assim, um processo y[n] passa a ter um modeloMA se for descrito pela Equacao 3.21, possuındo espectro de potencia e funcao densi-dade espectral de potencia apresentados nas Equacao 3.22 e 3.23 respectivamente. AFigura 3.3 mostra um diagrama do modelo MA.

y[n] =

q∑k=0

b[k] x[n− k] (3.21)

H(z) =

q∑k=0

b[k] z−k (3.22)

Py(ejw) = σ2

x

∣∣∣∣∣q∑

k=0

b[k] e−jwk

∣∣∣∣∣2

(3.23)

Os processos MA sao tambem conhecidos por modelos tudo-zero (all-zero). Ja osmodelos autorregressivos, tambem chamados tudo-polo (all-poles), sao casos especiaisdo ARMA onde os coeficientes b[k] sao iguais a zero, menos b[0] = 1, como se podever na Equacao 3.24. O espectro de potencia e a funcao densidade de potencia saodefinidas nas Equacoes 3.25 e 3.26 respectivamente. Como se pode ver, o PSD de umprocesso AR de ordem p possui 2p polos e nenhum zero, a nao ser aqueles localizadosem z = 0 e z = ∞. A Figura 3.4 ilustra este tipo de modelo.

y[n] = −p∑

k=1

a[k] y[n− k] + x[n] (3.24)

H(z) =1

1 +

p∑k=1

a[k] z−k

(3.25)

Py(ejw) =

σ2x∣∣∣∣∣1 +

p∑k=1

a[k] e−jwk

∣∣∣∣∣2 (3.26)

32


Figura 3.3: Diagrama de um modelo media movel (MA) para um processo estocastico,considerando a entrada x[n] e a saıda y[n].

Figura 3.4: Diagrama de um modelo autorregressivo (AR) para um processo es-tocastico, considerando a entrada x[n] e a saıda y[n].

33


Como ilustracao, apresenta-se um exemplo utilizando um filtro AR de segundaordem, que possui dois polos complexos conjugados, z1 = rejwo e z2 = re−jwo , sendowo igual a π/2, como apresentado pela Equacao 3.27. Assim, com base na Equacao3.25, encontra-se a funcao de densidade espectral de potencia para esse exemplo,considerando-se a variancia do ruıdo de entrada, σ2

x = 1, Equacao 3.28.

H(z) =1

(1− z1z−1)(1− z2z−1)=

1

1− (z1 + z2)z−1 + (z1z2)z−2(3.27)

Py(ejw) =

1

1 + r4 + 2r2cos(2w)(3.28)

A Figura 3.5 apresenta as densidades espectrais de potencia para tres casos, ondeos raios dos polos sao 0, 80; 0, 85 e 0, 90. Como se pode ver, conforme o raio vai seaproximando do cırculo unitario, os picos no espectro aumentam de tamanho e a sualargura de banda diminui.

Para o caso do modelo MA, a funcao PSD e o inverso da apresentada na Equacao3.28. Assim, conforme os modulos dos zeros se aproximam do cırculo unitario, maisprofundo sera o vale no espectro. A Figura 3.6 apresenta a funcao PSD para tres casos,com raios iguais a 0,75; 0,85 e 0,95.

Pode-se relacionar os modelos AR, ARMA e MA utilizando-se o teorema daDecomposicao WOLD, o qual diz que qualquer processo estocastico WSS pode serdecomposto em duas componentes, uma totalmente aleatoria e outra determinıstica.Um processo determinıstico e aquele que pode ser perfeitamente predito utilizando-seinfinitos dados passados.

Um exemplo desse teorema e um processo composto por ruıdo branco e umasenoide com fase aleatoria, sendo esta ultima caracterıstica que define o processo comoWSS [19]. Nesse caso, a parte aleatoria seria o ruıdo branco e a determinıstica asenoide. De forma semelhante, pode-se dizer que isso e equivalente a separar a funcaode densidade espectral de potencia, PSD (Power Spectral Density), em duas, uma partecontınua e outra discreta, devido ao ruıdo e a senoide respectivamente.

Se o PSD for puramente contınuo, qualquer processo AR ou ARMA pode serrepresentado por um unico modelo MA de ordem infinita. Dessa forma, caso o modeloescolhido para descrever um processo esteja errado, ainda se pode obter uma boaaproximacao usando uma ordem elevada.

Equacoes de Yule-Walker

As equacoes de Yule-Walker sao relacoes entre os parametros do modelo usadopara descrever um processo estocastico e a sua funcao de autocorrelacao. Para seencontrar essas relacoes, parte-se da Equacao 3.16 do modelo ARMA, multiplicando-seos dois lados por y∗[n − s] e tirando-se o valor esperado, obtendo-se a expressao aseguir, lembrando que ry[s] = Ey[n] y∗[n− s]:

ry[s] = −p∑

k=1

a[k] ry[s− k] +

q∑k=0

b[k] Ex[n− k] y∗[n− s]

34


Figura 3.5: Efeito da variacao do raio dos polos complexos conjugados de um filtroAR de segunda ordem, percebendo-se o aumento do pico e a diminuicao da largura debanda proporcionalmente a aproximacao do raio ao cırculo unitario.

Figura 3.6: Efeito da variacao do raio dos zeros complexos conjugados de um filtro MAde segunda ordem, percebendo-se o aumento da profundidade do vale e a diminuicaoda largura de banda proporcionalmente a aproximacao do raio ao cırculo unitario.

35


Como x[n] e WSS, se σ2y < ∞, y[n] tambem sera WSS [13], ou seja, considera-se

que o filtro seja estavel. Entao, x[n] e y[n] sao considerados processos conjuntamenteWSS (jointly WSS ). Assim, a esperanca do lado direito da expressao anterior e rxy[s−k]. Escrevendo-se essa correlacao cruzada em termos da autocorrelacao, ry[s], e daresposta impulsiva do filtro, h[n], tem-se a relacao a seguir:

rxy[s] = σ2x h∗[k − s]

Substituındo-se essa relacao na expressao anterior, chega-se as equacoes de Yule-Walker para o modelo ARMA, Equacao 3.29, sendo b[n] a resposta impulsiva de B(z).

ry[s] +

p∑k=1

a[k] ry[s− k] =

σ2x

q−s∑k=0

b[k + s] h∗[k], 0 ≤ s ≤ q

0, s > q

(3.29)

Com base nessas equacoes, tem-se uma relacao entre os coeficientes do filtro ea funcao de autocorrelacao, apesar desta ser nao-linear para os coeficientes devido aoproduto das respostas impulsivas h[n] e b[n].

Para o modelo AR, faz-se a substituicao: b[k] = δ[k]. Dessa forma, o somatoriodo lado direito de 3.29, modificado para 0 ≤ s ≤ q, resulta em h∗[−s], mas como essetermo e igual a zero para s > 0 e 1 para s = 0, chega-se as equacoes de Yule-Walkerdesse modelo, Equacao 3.30.

ry[s] +

p∑k=1

a[k] ry[s− k] = σ2x δ(s), s ≥ 0 (3.30)

Entao, dada a funcao de autocorrelacao, pode-se determinar os parametros domodelo AR, resolvendo-se um conjunto de equacoes lineares. Para o modelo MA,procede-se da mesma forma, mas as relacoes nao sao lineares. Para melhor visualizacaodesse resultado, apresenta-se as equacoes de Yule-Walker na forma matricial, Equacao3.31.

~Ry · ~a = ~ry, 1 ≤ s ≤ p (3.31)

ry[0] ry[−1] . . . ry[−p]ry[1] ry[0] . . . ry[−p + 1]

......

. . ....

ry[p] ry[p− 1] . . . ry[0]

1a[1]...

a[p]

=

σ2

x

0...0

A escolha do tipo de modelo a ser usado para descrever um processo e de grande

importancia, entao deve-se ter algumas ideias sobre o comportamento de cada um.O modelo MA possui PSD com banda larga, apresentando nulos quando os zeros dafuncao de sistema tem modulos proximos do cırculo unitario. Entao, utiliza-se essemodelo quando se precisa de vales no espectro. Ao contrario, quando o espectro possuipicos, deve-se usar um modelo AR, pois este possui polos, que levam a resposta parainfinito, tendo largura de banda estreita para ordens altas. Quando se precisam devales e picos, utiliza-se o modelo ARMA.

36


Quando se usa um modelo com caracterısticas diferentes daquelas do processoa ser descrito, os resultados serao pobres, mas podendo ser melhorados conforme seaumente a ordem. Por exemplo, caso se use um modelo MA para descrever um processocom PSD composto por picos, a falta de polos do modelo deve ser suprida com umamaior ordem.

Como ilustracao, apresentam-se dois exemplos utilizando um processo AR desegunda ordem sendo aproximados por modelos MA de ordens 2, 5 e 10. No primeirocaso, o modelo AR possui PSD dado pela Equacao 3.28, sendo o raio igual a 0,7071,que possui banda mais larga; e o segundo, com raio igual a 0,9592, de banda estreita.As aproximacoes pelos modelos MA sao apresentados para esses dois casos nas Figuras3.7 e 3.8 respectivamente.

Como se pode observar, no primeiro caso, o modelo MA de ordem 2 e 5 naoaproximam bem o AR(2), mas o MA(10) ja e um bom modelo. No segundo caso,nenhum dos tres apresenta bons resultados, pois o pico apresenta amplitude bem alta,entao, somente uma ordem elevada conseguiria aproxima-lo.

Figura 3.7: Comparacao do PSD do modelamento de um processo AR de segundaordem com banda larga (raio dos polos complexos conjugados iguais a 0,7071) pormodelos MA de ordem 2, 5 e 10.

Uma consideracao importante e feita quando se quer modelar um processo ARsomado com ruıdo branco, pois o resultado se comporta como ARMA. Isso pode servisto, utilizando-se um AR(2), de PSD PAR(ejw), somando com um ruıdo de varianciaσ2

v como exemplo, o PSD do processo resultante, Py(ejw), pode ser escrito como:

37


Figura 3.8: Comparacao do PSD do modelamento de um processo AR de segundaordem com banda estreita (raio dos polos complexos conjugados iguais a 0,9592) pormodelos MA de ordem 2, 5 e 10.

Py(ejw) = PAR(ejw) + σ2

v

=σ2

x

|A(ejw)|2+ σ2

v

=σ2

x + σ2v |A(ejw)|2

|A(ejw)|2

= PARMA(ejw) .

Entao, o resultado corresponde ao PSD de um processo ARMA(2,2). Dessa forma,quando a variancia σ2

v tender a zero, o processo sera bem descrito por um modelo AR.Quando esse valor aumenta, o filtro ARMA apresenta zeros com modulos perto docırculo unitario, ou seja, vales mais profundos no espectro de potencia, o que nao ebem representado por um AR. Mas se a ordem do modelo for aumentada, a aproximacaomelhora.

Isso pode ser observado nas Figuras 3.9 e 3.10 que apresentam dois processosARMA(2,2) sendo modelados por AR de segunda e quinta ordens [19]. Na primeirafigura, os zeros possuem modulos pequenos, 0,5, o que corresponderia a somar umruıdo de pequena variancia a um processo AR(2). Ja na segunda figura, os modulossao 0,9, correspondendo a variancias altas. Percebe-se que, no primeiro caso, o modeloAR de segunda ordem apresenta grandes erros, mas o de quinta aproxima muito bemo processo. Na Figura 3.10, nenhum dos dois aproxima bem, mas o de quinta ordemapresenta erros menores.

38


Figura 3.9: Aproximacao de um processo ARMA(2,2), que corresponde a um AR(2)mais ruıdo branco, por modelos AR de segunda e quinta ordens. Os zeros daqueleprocesso sao 0, 5 e0,6 π e 0, 5 e−0,6 π, com modulos longe do cırculo unitario, ou seja, avariancia do ruıdo e baixa.

Figura 3.10: Aproximacao de um processo ARMA(2,2), que corresponde a um AR(2)mais ruıdo branco, por modelos AR de segunda e quinta ordens. Os zeros daqueleprocesso sao 0, 9 e0,6 π e 0, 9 e−0,6 π, com modulos perto do cırculo unitario, ou seja, avariancia do ruıdo e alta.

39


3.3 Filtragem Adaptativa

Filtros de Wiener

Os filtros de Wiener sao uma classe de filtros otimizados, lineares e discretos notempo. Nestes, a saıda y[n] e uma estimativa de uma resposta desejada, d[n], dados umconjunto de amostras de entrada, x[n], de forma a minimizar o valor medio-quadraticodo erro, e[n] = d[n]− y[n], [14]. A Figura 3.11 apresenta um diagrama de blocos geraldos filtros Wiener.

Figura 3.11: Diagrama de blocos geral dos filtros Wiener, sendo x[n] a entrada, y[n] asaıda, d[n] a resposta desejada, e[n] o erro de estimativa e w os coeficientes do filtro.

A saıda do filtro e dada pela convolucao da serie temporal de entrada e a respostaimpulsiva do filtro, ambos com valores complexos, conforme a Equacao 3.32, sendo M aordem do filtro. Considera-se que a entrada e a resposta desejada sao series temporaisde processos estocasticos conjuntamente WSS, com media nula [14].

y[n] =M−1∑k=0

w[k]∗x[n− k] (3.32)

Para otimizar o filtro, busca-se minimizar o erro, chamado funcao custo, Equacao3.33, no sentido medio quadratico. Entao, para se encontrar o mınimo de ξ, utiliza-se acondicao de gradiente nulo em relacao aos coeficientes do filtro, ~∇ξ = 0. O resultado,chamado de princıpio da ortogonalidade, e apresentado na Equacao 3.34, apos algumasmanipulacoes matematicas [14], sendo eo[n] o erro de estimativa quando o filtro operaem condicao otima.

40


ξ = Ee[n] e∗[n] = E|e[n]|2 (3.33)

Ex[n− k] e∗o[n] = 0, k = 0, 1, . . . ,M − 1 (3.34)

Na condicao de otimizacao, a saıda do filtro pode ser considerada como umaestimativa otima da resposta desejada, assim y[n] = d[n]. Entao, pode-se reescrever oerro como apresentado a seguir, e quando se acham os valores medios quadraticos daexpressao, encontra-se o erro medio quadratico normalizado, Equacao 3.35.

eo[n] = d[n]− d[n]

ξmin = σ2d − σ2

d

ε = 1−σ2

d

σ2d

(3.35)

Nessas equacoes, os valores de σ2d e σ2

dsao as variancias da resposta desejada e

estimada, pois o valor medio quadratico de um processo WSS de media nula e iguala sua variancia. Avaliando-se o ε, percebe-se que o valor esta entre 0 e 1, sendo que,quanto menor, melhor e a estimativa da resposta desejada.

Relacoes importantes para o filtro de Wiener sao as Equacoes de Wiener-Hopf,as quais relacionam a autocorrelacao da entrada e a correlacao cruzada da entradae resposta desejada. Essas equacoes sao desenvolvidas a partir de 3.32 e 3.34, comomostram as expressoes a seguir, sendo wo[i] os coeficientes otimos do filtro:

E

x[n− k]

(d∗[n]−

M−1∑i=0

wo[i] x∗[n− i]

)= 0

M−1∑i=0

wo[i] Ex[n− k] x∗[n− i] = Ex[n− k] d∗[n]

M−1∑i=0

wo[i] r[i− k] = p[−k] (3.36)

Sendo r[i−k] e p[−k], a autocorrelacao de x[n] para um atraso k−i e a correlacaocruzada de x[n] e d[n] para um atraso de k respectivamente. Considerando o processoestocastico de entrada, ~x[n], e nao somente a serie temporal, x[n], encontra-se a formamatricial das Equacoes de Wiener-Hopf, 3.37, como apresentado na sequencia a seguir:

~x[n] = x[n], x[n− 1], . . . , x[n−M + 1]T

~R = E~x[n] ~xH [n]~p = E~x[n] d∗[n] = p(0), p(−1), . . . , p(1−M)T

~wo = wo[0], wo[1], . . . , wo[M − 1]T

~R ~wo = ~p ⇒ ~wo = ~R−1 ~p (3.37)

41


Predicao Linear

A predicao linear consiste em usar amostras de um processo estocastico WSS parafazer uma estimativa de outra amostra. Divide-se em dois tipos: predicao para frente(forward prediction) e para tras (backward prediction). No primeiro, usam-se M amos-tras passadas para prever a proxima, e no segundo, usam-se as M amostras presentespara prever a que vem antes delas. As Figuras 3.12 e 3.13 apresentam diagramas deblocos desses dois tipos de preditores respectivamente.

Figura 3.12: Diagrama de blocos de um preditor linear para frente, sendo x[n −1], . . . , x[n−M ] as entradas e x[n] a estimativa de x[n].

Com base na Figura 3.12, percebe-se que o preditor linear para frente e compostopor M entradas amostradas de um processo estocastico WSS de media nula, x[n −1], . . . , x[n − M ], que sao filtradas pelos coeficientes w∗

f [1], . . . , w∗f [M ], os quais sao

otimizados de forma media quadratica. Dessa forma, o valor estimado de x[n] e dadopela Equacao 3.38.

x[n] =M∑

k=1

w∗f [k] x[n− k] (3.38)

O erro de predicao e definido da mesma forma que para o filtro de Wiener,fM [n] = x[n]− x[n], e o valor medio-quadratico mınimo do erro e PM = E|fM [n]|2.Considerando o processo de entrada como ~x[n − 1] = x[n − 1], . . . , x[n − M ]T , oscoeficientes como ~wf = wf [1], . . . , wf [M ]T , a correlacao cruzada entre a entrada, aresposta desejada como ~r = E~x[n− 1] x∗[n] = r[−1], . . . , r[−M ]T e a variancia dex[n] como r[0], tem-se as Equacoes de Wiener-Hopf, Equacoes 3.39 e 3.40.

~R ~wf = ~r (3.39)

PM = r[0]− ~r H ~wf (3.40)

42


Figura 3.13: Diagrama de blocos de um preditor linear para tras, sendo x[n], . . . , x[n−M + 1] as entradas e x[n−M ] a estimativa de x[n−M ].

Quando o filtro preditor linear esta otimizado de forma media-quadratica, seuscoeficientes sao iguais aos do modelo autorregressivo. Combinando-se as duas equacoesanteriores, encontram-se as Equacoes de Wiener-Hopf aumentadas, Equacao 3.41, sendo~0, o vetor nulo com dimensoes M por 1.[

r[0] ~r H

~r ~R

] [1

−~wf

]=

[PM

~0

](3.41)

O preditor linear para tras e bem parecido com anterior, como se percebe naFigura 3.13. Este tambem e composto por M entradas amostradas de um processoestocastico WSS de media nula, x[n], . . . , x[n − M + 1], as quais sao filtradas peloscoeficientes w∗

b [1], . . . , w∗b [M ], otimizados de forma media quadratica. O objetivo e

estimar o valor de x[n − M ], conforme a Equacao 3.42, sendo o erro de estimativadefinido pela Equacao 3.43.

x[n−M ] =M∑

k=1

w∗b [k] x[n− k + 1] (3.42)

bM [n] = x[n−M ]− x[n−M ] e PM = E|bM [n]|2 (3.43)

Considerando ~wb = w∗b [1], . . . , w∗

b [M ]T , tem-se a correlacao cruzada entre aentrada e a resposta desejada como ~r B∗ = Ex[n] x∗[n − M ] = r(M), . . . , r[1]T ,sendo o operador .B correspondente ao arranjo ao contrario do vetor aplicado. Dessaforma, as Equacoes de Wiener-Hopf aumentadas sao definidas na Equacao 3.44. Arelacao entre os coeficientes dos preditores para frente e para tras e ~w B∗

b = ~wf .[~R ~r B∗

~r BT r[0]

] [−~wb

1

]=

[~0

PM

](3.44)

43


Pode-se definir os filtros de erro de predicao para frente e para tras a partir dasEquacoes 3.45 e 3.46, sendo os coeficientes aM [k] dados pela Equacao 3.47. Uma formade se calcular os coeficientes dos filtros de erro de predicao e usando o algoritmo deLevinson-Durbin [14], o qual utiliza a solucao das Equacoes de Wiener-Hopf aumenta-das de um filtro de ordem m− 1 para encontrar as solucoes de um de ordem m.

fM [n] =M∑

k=0

a∗M [k] x[n− k] (3.45)

bM [n] =M∑

k=0

a∗M [M − k] x[n− k] (3.46)

aM [k] =

1, k = 0

−wf [k], k = 1, 2, . . . ,M(3.47)

Sendo ~am−1 os coeficientes de um filtro de erro de predicao para frente de ordemm − 1, e ~aB∗

m−1 os coeficientes do filtro para tras de mesma ordem, pode-se encontraros novos coeficientes para um filtro de ordem m, utilizando-se as relacoes apresentadasnas Equacoes 3.48 e 3.49, respectivamente para os filtros para frente e para tras.

~am =

[~am−1

0

]+ κm

[0

~aB∗m−1

](3.48)

~aB∗m =

[0

~aB∗m−1

]+ κ∗m

[~am−1

0

](3.49)

As constantes κm sao chamadas coeficientes de reflexao. Partindo-se da Equacao3.48, faz-se a multiplicacao por ~Rm+1, que e a matriz de correlacao das entradas de umfiltro de ordem m. Assim, usando as Equacoes de Wiener-Hopf aumentadas, encontra-se uma relacao geral para o filtro para frente apresentada na Equacao 3.50, sendo ~0m,o vetor nulo de dimensoes m por 1.

[Pm

~0m

]=

Pm−1

~0m−1

∆m−1

+ κm

∆∗m−1

~0m−1

Pm−1

(3.50)

Dessa forma, obtem-se a definicao dos coeficientes de reflexao, Equacao 3.51, emfuncao do escalar ∆m−1, definido na Equacao 3.52. Encontra-se, tambem, a relacao derecorrencia para o erro de predicao mınimo, Equacao 3.53.

κm = −∆m−1

Pm−1

(3.51)

∆m−1 =m−1∑s=0

r[s−m] am−1[s] (3.52)

Pm = Pm−1 (1− |κm|2) (3.53)

A forma mais comum de se aplicar o algoritmo de Levinson-Durbin para encontraros coeficientes ~aM = aM [0], . . . , aM [M ] de um filtro de erro de predicao e utilizando a

44


funcao de autocorrelacao da entrada r[0], . . . , r[M ]. A partir desses valores, acham-se ∆m−1 e Pm−1 recursivamente, iniciando com os valores P0 = r[0] e ∆0 = r∗[1].Obtem-se os coeficientes de reflexao κm e, com eles, os coeficientes ~am. Termina-sequando m = M , sabendo que am[0] = 1 e am[k] = 0, para k > m. As estimativasdos coeficientes do filtro e do erro mınimo encontradas sao chamadas de estimativas deYule-Walker [14].

Filtros LMS

Os filtros de Wiener consideram que os sinais sao estacionarios (WSS), mas ge-ralmente esse nao o caso, assim outras tecnicas devem ser usadas. O metodo maiscomum para esse tipo de aplicacao e a filtragem LMS (Least-Mean-Square), perten-cente a famılia dos algoritmos adaptativos recursivos de gradiente estocastico [14]. Ofiltro LMS e bem simples, pois nao usa funcoes de correlacao ou inversoes matriciais.

A derivacao do algoritmo comeca pela aplicacao do contexto nao-estacionario aofiltro de Wiener. Dessa forma, considera-se que x[n] e d[n] nao sao conjuntamenteWSS na Equacao 3.32, entao a solucao nao e mais as Equacoes de Wiener-Hopf. Nocaso nao-estacionario, os coeficientes que minimizam o erro medio quadratico seraodependentes do tempo n, assim, tem-se a Equacao 3.54 para a estimativa da respostadesejada.

d[n] =M−1∑k=0

wn[k]∗x[n− k] = ~w Tn ~x[n] (3.54)

O projeto de um filtro adaptativo e bem mais difıcil que o de um invariante notempo, pois seria necessario achar um conjunto de coeficientes otimos para cada valorde n. Mas, isso e inviavel, assim, utiliza-se um fator de correcao aplicado a estimativados coeficientes para cada novo n [13], ∆~wn, conforme a Equacao 3.55. O diagrama deblocos de filtros LMS e igual ao do filtro de Wiener.

~wn+1 = ~wn + ∆~wn (3.55)

A regra usada para definir o fator de correcao e que o erro medio quadraticodeve diminuir com o tempo. Uma condicao e que, em um ambiente estacionario, oscoeficientes devem convergir para a solucao das equacoes de Wiener-Hopf. Entao, comessas consideracoes, utiliza-se a condicao de gradiente do erro igual a zero, encontrando-se as novas equacoes de Wiener-Hopf, Equacao 3.56, sendo ~R[n] uma matriz hermitianade autocorrelacoes de x[n] com dimensao M por M, e ~p[n] um vetor de correlacoescruzadas de x[n] e d[n] com dimensao M.

~R[n]~wn = ~p[n] (3.56)

Um caminho para a estimativa dos coeficientes e conhecido como metodo doGradiente Descendente, no qual o fator de correcao e a multiplicacao de um valorpositivo constante chamado de passo de estimativa, µ, e o negativo do gradiente do

45


erro, como apresentado na Equacao 3.57. Com alguma manipulacao matematica [13],encontra-se uma nova relacao, Equacao 3.58.

~wn+1 = ~wn − µ ~∇ξ[n] (3.57)

~wn+1 = ~wn − µ Ee[n] ~x ∗[n] (3.58)

O problema desse metodo e que a esperanca que aparece na expressao anteriore geralmente desconhecida, tendo que ser estimada, como, por exemplo, pela mediade amostra, Equacao 3.59. Entao, a partir desse ponto, introduz-se o algoritmo LMSque aproxima o somatorio somente pelo primeiro termo N = 1, resultando na Equacao3.60.

~wn+1 = ~wn −µ

N

N−1∑k=0

e[n− k] ~x ∗[n− k] (3.59)

~wn+1 = ~wn − µ e[n] ~x ∗[n] (3.60)

Com essa estimativa do gradiente, os coeficientes convergem para a solucao dasequacoes de Wiener-Hopf na media, pois o estimador e nao-influenciado (unbiased),

[13], conforme mostra a Equacao 3.61, sendo ~∇ξ[n] = −e[n] ~x ∗[n], o estimador.

E ~∇ξ[n] = −Ee[n] ~x ∗[n] = ~∇ξ[n] (3.61)

Para que o algoritmo LMS convirja na media, e necessario impor algumas res-tricoes no passo de estimativa. Inicialmente, considera-se 0 < µ < 2/λmax, sendo λmax

o maior autovalor da matriz de autocorrelacao de x[n], mas essa condicao nao poerestricoes na variancia dos coeficientes, assim, nao garante estabilidade em todos oscasos. Alem disso, e necessario se calcular a matriz de autocorrelacao.

Para melhorar essa condicao, utiliza-se o traco da matriz de autocorrelacao, tr(~R),como limite superior do maior autovalor. Com a consideracao que x[n] e WSS, tem-se

que tr(~R) = M r[0] = M E|x[n]|2, ficando a condicao conforme a Equacao 3.62.Nessa, utiliza-se E|x[n]|2 que representa a potencia de x[n], a qual e mais facil deser calculada.

0 < µ <2

M E|x[n]|2(3.62)

Uma medida de convergencia do algoritmo LMS e o desajuste (adjustment), M,o qual e considerado como um erro medio-quadratico normalizado, definido aproxima-damente, pela Equacao 3.63, quando o passo e µ << 2/λmax.

M≈ 1

2µ tr(~R) (3.63)

O algoritmo LMS e muito usado para se estimar os coeficientes de preditoreslineares. Como ilustracao, apresenta-se um exemplo de processo autorregressivo, geradopela Equacao 3.64, sendo x[n] um processo de ruıdo branco de media nula e varianciaunitaria. A funcao de sistema e dada pela Equacao 3.65. Para se encontrar o preditor

46


linear otimo para esse processo, Equacao 3.66, utiliza-se o algoritmo LMS de segundaordem. Os valores esperados para os coeficientes sao 0, 75 e −0, 9.

y[n] = 0.75y[n− 1]− 0.9y[n− 2] + x[n] (3.64)

H(z) =1

1− 0.75z−1 + 0.9z−2(3.65)

y[n + 1] = wn[1] y[n] + wn[2] y[n− 1] (3.66)

O sinal desejado e encontrado deslocando-se y[n] em uma amostra para a frente,assim d[n] = y[n + 1]. O valor do passo e encontrado usando a Equacao 3.62, queresultou em µ = 0.195, mas para se garantir a convergencia, usaram-se valores bemmenores: 50 vezes menor para um primeiro caso e 100 vezes para outro. Quanto menoro valor do passo, menor sera o desajuste mas, mais lenta a conversao. A Figura 3.14apresenta os resultados das simulacoes dos dois casos, sendo o primeiro apresentadoem cima e o segundo embaixo.

Figura 3.14: Valores dos coeficientes para o modelo LMS aplicado a um processoautorregressivo de segunda ordem. O grafico de cima apresenta os resultados parapasso dividido por 50 e, o de baixo, para passo dividido por 100.

47


Para se achar o valor dos coeficientes, encontrou-se a media dos 256 ultimos pontosdas simulacoes, pois, como se pode ver, estes flutuavam ao redor do valor esperado(linhas pontilhadas). No primeiro caso, encontrou-se w[1] = 0, 786 e w[2] = −0, 862,para o segundo, w[1] = 0, 778 e w[2] = −0, 874.

Metodo de Prony e da Autocorrelacao

Uma tecnica muito usada para se estimar os coeficientes de um modelo AR apartir de um processo y[n] de tamanho finito N e o metodo da autocorrelacao. Essatecnica se baseia no metodo de Prony, o qual busca a solucao de um conjunto deequacoes lineares, como o algoritmo de Levinson-Durbin.

No metodo de Prony, um sinal x[n] e modelado, de forma geral, como a respostaimpulsiva de um filtro linear invariante no tempo com funcao de sistema igual a ummodelo ARMA. Com a consideracao que x[n] = 0 para n < 0, busca-se os coeficientesa[k] e b[k] que minimizam o erro quadratico E , definido a partir do erro e′[n], Equacao3.67, sendo h[n] a resposta impulsiva do filtro como mostra a Figura 3.15. A Equacao3.68 apresenta a resposta em frequencia do erro, E ′(z).

e′[n] = x[n]− h[n] (3.67)

E ′(z) = X(z)− B(z)

A(z)(3.68)

Figura 3.15: Diagrama de blocos do metodo de prony, mostrando a forma de se encon-trar o erro e[n], sendo b[n] a estimativa de b[n].

Multiplicando-se os dois lados da Equacao 3.68 por A(z), encontra-se um novoerro, definido, em z, por E(z) = A(z) E ′(z), assim, o erro no domınio do tempo eapresentado pela Equacao 3.69, sendo a[n] e b[n] as respostas impulsivas de A(z) eB(z) e b[n] = 0 para n > q. Assim, o erro quadratico fica como na Equacao 3.70.

e[n] = a[n] ∗ x[n]− b[n] (3.69)

E =∞∑

n=q+1

|e[n]|2 (3.70)

Considerando o princıpio da ortogonalidade apresentado para Filtros de Wienerna Equacao 3.34, e o erro e[n] da Equacao 3.67, chega-se a Equacao 3.71, para n > q,

48


quando b[n] = 0. Com algumas substituicoes e considerando-se a autocorrelacao r[k, l]apresentada na Equacao 3.72 [13], chega-se nas Equacoes 3.73, chamadas de Equacoes

Normais de Prony, sendo ~R, a matriz composta pelas autocorrelacoes r[1, 1] a r[p, p] e~r = r[1, 0], . . . , r[p, 0]T .

∞∑n=q+1

e[n] x∗[n− k] = 0, k = 1, 2, . . . , p (3.71)

r[k, l] =∞∑

n=q+1

x[n− l] x∗[n− k] (3.72)

Percebe-se que essas equacoes possuem a mesma forma das Equacoes de Wiener-Hopf, mas usando o vetor ~a = a[1], a[2], . . . , a[p], semelhante a definicao da Equacao3.47 da Predicao Linear.

p∑l=1

a[l] r[k, l] = −r[k, 0], k = 1, 2, . . . , p (3.73)

~R~a = −~r

Considerando-se a matriz ~Xq composta por p colunas infinitas, como apresentadona Equacao 3.74, pode-se reescrever as Equacoes Normais como na Equacao 3.75, sendo~xq+1 = x[q + 1], x[q + 2], . . .T .

~Xq =

x[q] x[q − 1] . . . x[q − p + 1]

x[q + 1] x[q] . . . x[q − p + 2]x[q + 2] x[q + 1] . . . x[q − p + 3]

......

...

(3.74)

( ~X Hq

~Xq) ~a = − ~X Hq ~xq+1 (3.75)

Quando se considera um modelo composto somente por polos (all-pole), como oAR, pode-se fazer algumas simplificacoes no metodo de Prony. Inicialmente, modifica-se os limites do somatorio do erro quadratico, 3.70, pois q = 0, incorporando-se o casode n = 0 tambem [13], Equacao 3.76. O erro e[n] pode ser simplificado a partir daEquacao 3.67, retirando b[n] e lembrando que a[0] = 1, resultando na Equacao 3.77.

E =∞∑

n=0

|e[n]|2 (3.76)

e[n] = x[n] +

p∑k=1

a[k]x[n− k] (3.77)

Considerando-se que a autocorrelacao definida na Equacao 3.72 tem o limiteinferior do somatorio modificado para n = 0, pode-se tambem mudar a notacao der[k, l] para r[k− l], pois, no caso de filtros all-pole, essa funcao so depende da diferenca

49


dos coeficientes, como em processos WSS [13]. Assim, tem-se as Equacoes Normaisdesses filtros, Equacao 3.78.

p∑l=1

a[l] r[k − l] = −r[−k], k = 1, 2, . . . , p (3.78)

O metodo da Autocorrelacao e baseado no Prony, sendo utilizado para casos ondeo processo de entrada x[n] e finito, conhecido no intervalo [0, N ]. Entao, aplica-se umajanela temporal a x[n], de forma a o tornar nulo fora do intervalo conhecido, para, emseguida, achar os coeficientes do modelo.

Para isso, busca-se minimizar o erro quadratico, E , da mesma forma que nometodo de Prony. Como x[n] so e conhecido no intervalo [0, N ], nao se pode encontraro erro para n < p e n > N , entao, faz-se x[n] = 0 para n < 0 e n > N , criando umnovo sinal, x[n], pela aplicacao de uma janela retangular, wR[n], como apresentado naEquacao 3.79. Assim, para o calculo da autocorrelacao de x[n], utiliza-se a Equacao3.80.

x[n] = x[n] wR[n] (3.79)

wR[n] =

1, n = 0, 1, . . . , N0, caso contrario

r[k] =∞∑

n=0

x[n] x ∗[n− k] =N∑

n=k

x[n] x∗[n− k], k = 0, 1, . . . , p (3.80)

As Equacoes Normais para o caso de modelos all-pole do metodo Prony se apli-cam a este metodo, podendo-se resolve-las usando o algoritmo de Levinson-Durbinapresentado. Apesar do uso da janela acrescentar um erro a estimativa dos coeficien-tes, tem-se a vantagem de garantir a estabilidade do modelo, pois os polos de H(z)estarao dentro do cırculo unitario [13].

Pode-se tambem expressar o metodo da autocorrelacao como uma busca dasolucao de mınimos quadraticos (least squares) de um conjunto de equacoes linea-res sobre-determinadas. Procura-se por um conjunto de coeficientes a[k], de forma

que o erro e[n] = 0 para n > 0, Equacao 3.81, sendo ~Xp uma matriz com di-mensoes (N + p) por p, formada por versoes deslocadas da coluna x[0], x[1], . . . , x[p−1], x[p], . . . , x[N − 2], x[N − 1], x[N ], 0, . . . , 0T e ~x1 = x[1], x[2], . . . , x[N ], 0, . . . , 0T ,com dimensao (N + p). A solucao de mınimos quadraticos e determinada pela solucaodas Equacoes lineares 3.82.

p∑k=1

a[k]x[n− k] = −x[n] ⇒ ~Xp ~a = −~x1 (3.81)

( ~X Hp

~Xp) ~a = ~Xp ~x1 (3.82)

50

3.4 Rastreamento de Polos


Quando se consideram sinais nao-estacionarios, muitos metodos para obtencao do es-pectro de potencia se tornam incapazes de descrever informacoes importantes do sinal,entao, devem-se usar metodos mais refinados para a estimativa do PSD [24, 25]. Nestasecao, apresenta-se uma tecnica que usa o rastreamento de polos da decomposicaoautorregressiva do sinal para estimar caracterısticas interessantes do espectro.

O modelo AR apresenta a vantagem da decomposicao espectral rapida pelo usode um algoritmo residual [24]. Utilizando estimativa dos polos, divide-se o espectro emcurvas com formato de sino, sendo as caracterısticas de potencia e frequencia de cadacomponente espectral obtidos pela posicao e resıduo de cada polo.

Apesar dessa tecnica apresentar bons resultados e ser bem estabelecida, ela naoincorpora a extracao e rastreamento dos polos, os quais devem ser obtidos dos coefi-cientes do modelo. Segundo Mainardi [24], os metodos tradicionais de fatoracao paraestimativa dos polos sao ineficientes para tal, porque precisam ser reiniciados cada vezque novos coeficientes sao apresentados, nao usando os valores dos polos estimadosanteriormente.

Dessa forma, apresentam-se dois metodos de rastreamento dos polos que em-pregam as estimativas anteriores para atualizar os valores a cada iteracao, estes sao:metodo da linearizacao e de Bairstow. No Capıtulo 4, introduz-se outro metodo ba-seado nos autovalores da matriz companheira dos coeficientes, que apesar de nao usaras estimativas anteriores, apresentou bom desempenho nas simulacoes. Inicialmente,apresentam-se os conceitos que relacionam os polos da decomposicao AR com o espectrodo sinal.

Alguns trabalhos da area de bioengenharia apresentaram bons resultados como uso dos polos de modelos autorregressivos para caracterizacao de sinais. Pode-secitar, Patomaki [32], que rastreou um polo de um modelo ARMA(4,2) de sinais nao-estacionarios de EEG (Eletroencefalograma) para detectar picos correspondentes a mu-dancas de estados de sono em ratos. Para se obter o polo, usou-se uma iteracao dometodo de Newton.

Em [25], utiliza-se o rastreamento de polos do modelo AR para se encontrar oespectro de tempo-frequencia em tempo real usando o Metodo de Bairstow. Com isso,foi possıvel obter mais informacoes sobre o status do sistema nervoso autonomo depacientes com doenca na arteria coronaria que sofriam isquemia induzida por drogas.

Em [6], procura-se estimar o espectro de potencia de um sinal de voz, o quale geralmente modelado por linhas espectrais discretas correspondendo as harmonicasde pitch. Assim, utiliza-se o modelo AR para essa estimativa, pois este apresenta oenvelope do espectro mais proximo do esperado para os sinais usados.

Em [28], usa-se um modelo ARMA para estimar o espectro de potencia de sinaisde ECG (Eletrocardiograma), assim como em [24, 25], os quais usam modelos AR. Como espectro, foi possıvel complementar exames clınicos, levando a melhores diagnosticos.Porque o metodo apresenta qualidades mais interessantes que o FFT (Fast FourierTransform), o qual necessita de janelamento.

Conforme apresentado anteriormente, a Equacao 3.25 mostra o espectro de potenciaH(z) de um modelo AR, o qual se compoe de uma fracao com numerador unitario edenominador polinomial. Decompondo-se esse polinomio em um produtorio, chega-se

51


a uma forma alternativa do espectro de potencia, Equacao 3.83, sendo zk os polos deH(z).

H(z) =1

1 +

p∑k=1

a[k] z−k

=1

p∏k=1

(1− zk z−1)

(3.83)

Quanto a extracao dos componentes oscilatorios do sinal, a representacao porpolos oferece uma compreensao mais rapida dos fenomenos, devido a relacao unicaentre os polos e os picos espectrais [24]. Cada polo real e cada par de polos complexoscontribui para um pico em formato de sino no PSD. A frequencia central, fk, de cadapico e obtida da fase dos polos, φk, e da frequencia de amostragem do sinal, fs, conformea Equacao 3.84.

fk =φk

2πfs (3.84)

Uma caracterıstica muito importante dessa analise provem do fato de que quantomais proximo um polo estiver do cırculo unitario no plano z, maior sera o seu picocorrespondente no espectro [47]. Para melhor compreensao desses conceitos, apresenta-se um exemplo.

Neste exemplo, utiliza-se um processo AR dado pelas Equacoes 3.85, sendo x[n]um ruıdo branco de media nula e variancia σ2

x. Os polos de H(z) sao 0, 3626; −0, 8325+0, 4425i; −0, 8325−0, 4425i; 0, 5012+0, 6077i e 0, 5012−0, 6077i, apresentados no planoz na Figura 3.16 e o PSD de y[n] na Figura 3.17.

y[n]+0, 3x[n−1]−0, 4x[n−2]+0, 2x[n−3]+0, 5x[n−4]−0, 2x[n−5] = x[n] (3.85)

H(z) =1

1 + 0, 3z−1 − 0, 4z−2 + 0, 2z−3 + 0, 5z−4 − 0, 2z−5

Como se pode ver no espectro, Figura 3.17, tem-se dois picos mais um valornao-nulo de potencia na origem, os quais correspondem aos pares de polos complexose ao polo real. O par complexo do semiplano esquerdo do plano z possui frequencianormalizada igual a 0,4222 Hz, o outro par, 0,1402 Hz, que correspondem as frequenciascentrais dos picos. Como se pode ver, o par de polos que possui o modulo mais proximodo cırculo unitario, 0,9428, contribui com o maior pico no espectro.

Metodo da Linearizacao e de Bairstow

Os dois metodos apresentados a seguir sao: metodo da linearizacao, o qual utilizaum procedimento de linearizacao da relacao nao-linear entre os polos e os coeficientesdo modelo AR, apresentado em [24, 30], e metodo de Bairstow, o qual e um algoritmorecursivo muito usado para se calcular as raızes complexas de polinomios [24, 26]. Emambos, o procedimento e o mesmo, quando se obtem novos coeficientes do modelo,estes sao utilizados para estimar recursivamente as novas posicoes dos polos.

52


Figura 3.16: Polos de H(z) do exemplo apresentados no plano z, sendo 0, 3626;−0, 8325 + 0, 4425i; −0, 8325− 0, 4425i; 0, 5012 + 0, 6077i e 0, 5012− 0, 6077i.

Figura 3.17: Densidade espectral de potencia de y[n] mostrando dois picos, o menorcom frequencia central normalizada 0,1402 Hz, o maior com frequencia 0,4222 Hz.

53


No metodo da linearizacao, considera-se a Equacao 3.86, que relaciona os polose coeficientes do modelo AR. Para se verificar como sao afetadas as posicoes dos polosquando os coeficientes sofrem variacoes, recorre-se a um procedimento de linearizacao,sendo que a mudanca na posicao do i-esimo polo e definida conforme a Equacao 3.87,sendo a derivada δzi/δak definida na Equacao 3.88.

1

1 +

p∑k=1

a[k] z−k

=1

p∏i=1

(1− zi z−1)

(3.86)

∆zi =

p∑k=1

δzi

δak

∆ak (3.87)

δzi

δak

=zp−k

ip∏

l=1|l 6=i

(zi − zl)

(3.88)

Assim, encontrando-se o valor de ∆zi, pode-se estimar a nova posicao do polozi, sendo zi(t) = zi(t − 1) + ∆zi. Duas consideracoes devem ser feitas, primeira,essas equacoes se aplicam quando as mudancas nos coeficientes do modelo AR saopequenas, e segunda, o algoritmo e bastante sensıvel a grupos muito proximos depolos, apresentando erros consideraveis [24].

O metodo de Bairstow e um algoritmo que pode ser usado para encontrar raızescomplexas de polinomios utilizando somente aritmetica real, com a limitacao de os coe-ficientes serem reais, pois assim as raızes complexas ocorrem aos pares conjugados. Estemetodo baseia-se no algoritmo de deflacao quadratica. Considerando-se o polinomiop(x) de coeficientes reais e grau n, Equacao 3.89, pode-se encontrar um fator quadraticoq(x) = x2 − rx− s, reescrevendo-se o polinomio como p(x) = q(x) Q(x) + R(x), sendoQ(x) um polinomio de grau n − 2, Equacao 3.90, e R(x) um resto de grau menor ouigual a 1, Equacao 3.91.

p(x) = a1xn + a2x

n−1 + . . . + anx + an+1 (a1 6= 0) (3.89)

Q(x) = b1xn−2 + b2x

n−3 + . . . + bn−2x + bn−1 (3.90)

R(x) = bn(x− r) + bn+1 (3.91)

Expandindo a decomposicao de p(x), encontram-se as relacoes entre os coefici-entes an, bn, r e s, conforme a Equacao 3.92. Considerando o algoritmo da deflacaoquadratica, encontra-se o metodo de Bairstow, cujo objetivo e buscar os coeficientes re s de q(x) para os quais bn e bn+1 sejam iguais a zero, ou seja, o resto R(x) = 0. Dessaforma, q(x) apresenta duas raızes do polinomio p(x). Entao, pode-se retirar esse termodo polinomio e continuar a recursao em Q(x) para buscar as outras raızes.

bn = an + rbn−1 + sbn−2 (3.92)

Buscam-se os valores otimos r e s, partindo de valores iniciais r0 e s0, pela adicaode fatores corretivos, ∆rk e ∆sk, conforme as Equacoes 3.93, obtendo, a partir de

54

3.5 Simulacao de Monte Carlo

valores anteriores rk e sk, outros que estejam mais proximos do esperado, rk+1 e sk+1.

rk+1 = rk + ∆rk (3.93)

sk+1 = sk + ∆sk (3.94)

Para se encontrar esses fatores corretivos, faz-se uma nova consideracao, aposrealizar a divisao polinomial para encontrar os coeficientes bn, conforme 3.92, realiza-seuma nova divisao, de Q(x) por q(x), encontrando os coeficientes cn, que se relacionam abn da mesma forma que estes a an, (Equacao 3.92). Entao, conforme [26], encontram-seos fatores, Equacoes 3.95 e 3.96.

∆rk =bncn−1 − bn+1cn−2

cncn−2 − c2n−1

(3.95)

∆sk =bn+1cn−1 − bncn

cncn−2 − c2n−1

(3.96)

Os valores iniciais sao muito importantes para a convergencia rapida do metodo,normalmente, escolhe-se r0 e s0 iguais a zero. Mas, para se encontrar raızes grandes,pode-se usar r0 = −a2/a1 e s0 = −a3/a1, para pequenas, r0 = −an/an−1 e s0 =−an+1/an−1. Quando o q(x) estimado nao e totalmente exato, deve-se fazer uso dealgumas tecnicas para certificacao dos resultados [26].


O metodo de simulacao de Monte Carlo, segundo [36], baseia-se essencialmente nasimulacao de variaveis aleatorias. Sua origem se deve ao artigo The Monte CarloMethod publicado em 1949, de autoria de S. Ulam, com grandes contribuicoes de Johnvon Neumman.

O nome do metodo provem da cidade do principado de Monaco, chamada MonteCarlo, mundialmente conhecida por seus cassinos, os quais usam a roleta, um dosaparelhos mais simples para geracao de numeros aleatorios.

O algoritmo de calculo do metodo e simples e depende de cada aplicacao, consis-tindo basicamente em N simulacoes independentes da mesma experiencia e a obtencaode resultados estatısticos, variando somente os numeros aleatorios utilizados para asvariaveis. Por isso, a simulacao de Monte Carlo tambem e chamada de metodo dasprovas estatısticas [36].

Pode-se considerar que o erro inerente ao metodo e proporcional a√

D/N , sendoD uma constante dependente da aplicacao utilizada, a qual pode ainda ser modificadapara cada aplicacao para melhorar os resultados. Dessa forma, com este metodo nao seconsegue uma precisao elevada (acima de 95%) sem utilizar um numero muito grandede simulacoes, o que e muito custoso computacionalmente.

O metodo de Monte Carlo permite, entao, simular o comportamento de proces-sos que dependem de fatores aleatorios, contanto que se possa associar um modeloprobabilıstico ou um modelo artificial que aproximem esses processos [36].

55


Nas simulacoes, como ja mencionado, faz-se necessario suprir um conjunto denumeros aleatorios. Para isso, implementam-se algoritmos geradores de numeros pseu-doaleatorios, os quais sao obtidos a partir de uma prescricao (uma regra) e que possuemaltos perıodos sem repeticao, sendo dependentes de uma semente para serem iniciados.

Para ilustrar a utilizacao do metodo, apresenta-se um exemplo de estimativa dovalor de π. Considera-se uma circunferencia centrada na origem de raio unitario inscritaem um quadrado de lado dois. Geram-se pontos (x, y), sendo as coordenadas x e yduas variaveis aleatorias uniformemente distribuıdas entre 0 e 1. Levando-se em contasomente a regiao onde x e y tem valores entre 0 e 1, as areas das duas figuras sao 1e π/4, assim, sendo n o numero de pontos gerados e k o numero de pontos que caemdentro do cırculo, ou seja cujas coordenadas tenham a relacao

√(x2 + y2) ≤ 1, tem-se,

pela proporcao das areas, que k/n = π/4, assim, π = 4k/n.Para verificar o erro dessa estimativa, algumas simulacoes foram feitas utilizando-

se um gerador de numeros aleatorios padrao. A Figura 3.18 mostra um histograma dasdistancias dos pontos aleatorios ate a origem para o caso de n = 500.

Figura 3.18: Histograma das distancias dos pontos aleatorios ate a origem para o casode n = 500.

A Figura 3.19 apresenta o erro de estimativa para 100 simulacoes com diferentesvalores de n, comecando em 500 ate 50000. Percebe-se que o erro diminui com oaumento do numero de pontos.

Conforme o metodo de Monte Carlo, para se obter uma estimativa do resultado,no caso o valor estimado de π, deve-se fazer varias simulacoes. Dessa forma, para ocaso de n = 500, realizou-se N = 100 simulacoes, obtendo-se um erro de estimativamedio de 1, 81%, com desvio padrao de 1, 45%, o que esta de acordo com o esperadopelo metodo.

56


Figura 3.19: Erro de estimativa para 100 simulacoes com diferentes valores de n,comecando em 500 ate 50000. Percebe-se o erro diminui com o aumento do numero depontos.

57

58

Capıtulo 4

Materiais e Metodos

Neste capıtulo, serao apresentados os meios de teste usados, descrevendo a montageme utilizacao nos estudos de MSS. Em seguida, mostram-se os equipamentos e softwarede aquisicao de sinais de ultra-som.

Apos essas apresentacoes, serao descritos os softwares de simulacao de ecos deultra-som, explicando o gerador desses ecos, o qual e baseado nos modelos propostos;e de analise, mostrando todas as etapas envolvidas com exemplos para ilustracao.

Por fim, os aspectos computacionais dos algoritmos sao explicados, especifica-mente os metodos da convolucao de dois sinais, transformada de Hilbert, autocor-relacao, LMS e calculo de raızes complexas de polinomios usando estimativa de auto-valores da matriz companheira.

4.1 Phantoms

Para uma melhor compreensao dos processos de retroespalhamento das ondas de ultra-som em tecidos biologicos, utilizam-se meios de teste com propriedades semelhantes aesses tecidos. Esses meios sao chamados de phantoms, sendo construıdos segundo osmodelos propostos para os tecidos.

Os materiais mais usados para fabricacao dos phantoms sao alguns tipos de geis,gelatina ou agar, os quais simulam a consistencia de um tecido mole, e, para modelaros espalhadores, usa-se po de grafite ou cal [8]. Em [15], usou-se uma espuma feitade uma malha de poliuretano. Nesses trabalhos, o uso de meios de teste contribuiupara a compreensao dos modelos, que puderam ser melhorados para, em seguida, seremestudados em tecidos reais.

Os phantoms usados neste trabalho foram feitos com gelatina e po de grafite.Estes eram em formato cilındrico, feitos em tamanhos que variavam entre 1 cm e 4, 5 cmde altura, possuindo base circular de 3 cm de diametro. Para simular a parte difusa domodelo, utilizou-se o po de grafite, que ficou homogeneamente espalhado pelo meio.

Para simular o espacamento medio dos espalhadores, foram colocadas de duas aseis interfaces de tres tipos diferentes de materiais entre as camadas de gelatina. Essesmateriais eram: fios de cobre, folhas de acetato ou finas camadas feitas de gelatinacom o dobro da concentracao de grafite. Os espacamentos variavam de 0, 1 a 1 cm.Devido as suas constituicoes, essas interfaces gerariam ecos de maior amplitude quandofossem inspecionadas por ultra-som. A Figura 4.1 ilustra a constituicao dos phantoms

59

Capıtulo 4 - Materiais e Metodos

construıdos e a Figura 4.2 apresenta uma foto de tres phantoms usados para se coletarsinais, podendo-se ver as marcas de tres interfaces no phantom da esquerda e quatrono da direita, no do meio nao sao visıveis.

Figura 4.1: Ilustracao da constituicao dos phantoms usados como meios de teste parainspecao de ultra-som.

Figura 4.2: Foto de tres phantoms usados para se coletar sinais. Pode-se ver as marcasde tres interfaces no phantom da esquerda e quatro no da direita, no do meio nao epossıvel se observar.

4.2 Equipamentos e Software de Aquisicao

O equipamento usado no estudo consiste em um microcomputador Intel Pentium III,433 MHz, 128 megabytes de memoria RAM, com os sistemas operacionais MicrosoftWindows 98 e Redhat Linux 9.0, um sistema de aquisicao de dados e outro de posici-onamento.

O sistema de aquisicao de dados foi adquirido da empresa Optel, modelo OPKUD100, consiste em: uma placa com barramento ISA, que gera o sinal de excitacao dotransdutor e recebe o sinal de eco proveniente do mesmo; um modulo de potencia, que

60


amplifica o sinal de excitacao, com cabo serial para chegada do sinal da placa e cabocoaxial para transmissao do sinal de eco; dois transdutores, ambos feitos de ceramicapiezeletrica (PZT). Um e focalizado, com diametro de 10 mm e frequencia de operacaode 8 MHz; o outro, nao focalizado, com diametro de 15 mm e frequencia de 3, 5 MHz.A Figura 4.3 apresenta uma foto com os dois transdutores.

Figura 4.3: Os dois transdutores usados nas coletas de sinais de ultra-som. O dadireita e focalizado com diametro de 10 mm e frequencia de operacao de 8 MHz, e oda esquerda e nao focalizado com diametro de 15 mm e frequencia de 3, 5 MHz.

Figura 4.4: Ilustracao da montagem experimental usada nas coletas de sinais de ultra-som de meios de teste.

61


Figura 4.5: Foto do equipamento usado nos testes, montado no Laboratorio de Instru-mentacao Biomedica da UEL.

O sistema de posicionamento foi adquirido da empresa Servo Systems Co., con-siste em: uma placa de controle barramento ISA, modelo PRO500, ACSTech80, queenvia o sinal de posicao e recebe sinais de status do sistema; um modulo de potencia,modelo AMP2035, para amplificar o sinal e acionar o motor de passo; um eixo linear(posicionador) com motor de passo em rosca sem fim, modelo LPS-A4 e cabo de sinaisparalelo para conexao entre o modulo de potencia e o computador. Completando o sis-tema, usam-se uma fonte regulada de 5V e um suporte, usado para fixar o eixo linear(junto com o motor) e o transdutor.

Tambem fazem parte do equipamento: um tanque acustico, com dimensoes de 30por 15 cm de base e 20 cm de altura, feito de vidro e uma base de metal cilındrica, de 12cm de diametro e 1 cm de altura, feita de aco inox. A Figura 4.4 ilustra o equipamentousado nas coletas de dados e a Figura 4.5 apresenta uma foto do equipamento noLaboratorio de Instrumentacao Biomedica da UEL.

Alguns conceitos sobre o campo de ultra-som gerado por transdutores circularessao importantes na montagem experimental para se obter um bom sinal retroespalhadopela amostra. Sabe-se que o campo na frente do transdutor depende de seu raio,frequencia de ressonancia e do meio [11]. A onda de ultra-som entra no meio definindoduas regioes diferentes: a zona de Fresnel ou campo proximo, na qual ela e cilındrica,e a zona de Fraunhofer ou campo distante, longe do transdutor, onde a onda diverge

62


[20]. Na pratica, essa divisao nao e bem definida, existindo uma zona de transicao. AFigura 4.6 ilustra a disposicao dessas regioes [11].

No campo proximo, o ultra-som irradiado por diferentes partes do elemento pi-ezeletrico viaja como ondas esfericas que interferem construtiva e destrutivamente,portanto existem regioes de maximo e mınimo ao longo da onda. O comprimento L docampo proximo depende do raio r do transdutor e do comprimento de onda λ do somno meio na frente do transdutor, conforme a Equacao 4.1.

L =r2

λ(4.1)

Figura 4.6: Ilustracao da disposicao do campo de ultra-som criado por um transdutorcircular com impedancia casada com a do meio.

No campo distante, o ultra-som diverge e aparenta estar vindo de um pontolocalizado no centro do transdutor, como mostra a Figura 4.6. O angulo de divergenciaθ tambem depende do raio do transdutor e do comprimento de onda do som no meio,conforme a Equacao 4.2.

sen(θ) =0, 61λ

r(4.2)

Assim, o raio do transdutor e o comprimento de onda do ultra-som no meiodeterminam a natureza do campo em frente ao transdutor. O comprimento do campoproximo aumenta com o aumento do diametro ou frequencia de operacao do transdutor.A divergencia do feixe no campo distante decresce nessas mesmas condicoes. Essasdiscussoes se baseiam em transdutores com elementos piezeletricos planos, como o casodo transdutor de 3, 5 MHz utilizado nos experimentos.

Pode-se fazer um melhor controle do feixe de ultra-som utilizando-se transdutoresfocalizados, como o transdutor de 8 MHz utilizado nos experimentos. A focalizacaoe realizada com lentes acusticas, refletores curvos, ou mesmo, elemento piezeletricosde superfıcie curva. A Figura 4.7 apresenta os campos proximo e distante para umtransdutor focalizado no qual se utiliza lente [11], sendo o feixe convergente antes de

63


uma regiao de foco e divergente depois. Para o uso deste tipo de transdutor, deve-secolocar a amostra a ser inspecionada por ultra-som na regiao de foco.

Figura 4.7: Ilustracao dos campos proximo e distante para um transdutor focalizadono qual se utiliza lente, apresentando-se tambem a regiao de foco.

Entao, a regiao ideal para se colocar a amostra para o transdutor de 3, 5 MHz,que e plano, e no final do campo proximo, antes de se comecar o campo distante [9].Assim, considerando a velocidade do ultra-som media em tecidos moles como c = 1540m/s, sendo o raio igual a r = 7, 5 mm e o comprimento de onda λ = 0, 440 mm, tem-seuma distancia L = 12, 78 cm. O transdutor de 8 MHz e focalizado, sendo a distanciafocal igual a 10 cm.

O software utilizado para aquisicao de dados e chamado SADUS (Sistema deAquisicao de Dados de Ultra-Som), atualmente na versao 5.0, foi desenvolvido para usono Laboratorio de Instrumentacao Biomedica do Departamento de Engenharia Eletrica(DEEL) da Universidade Estadual de Londrina (UEL). Utilizaram-se as linguagensTcl/Tk e C++ para plataforma Redhat Linux, apresentando interface de usuario graficaGUI (Graphic User Interface) [42].

O SADUS faz a comunicacao do computador com os sistemas de aquisicao dedados e de posicionamento, usando o hardware de controle instalado com o proto-colo de comunicacao fornecido pelas empresas fabricantes. Com ele, pode-se ajustaros parametros da inspecao por ultra-som, apresentados na Tabela 4.1, os quais sao:tamanho do vetor de dados (buffer), ganho aplicado ao sinal recebido (gain), tensaoeletrica do pulso aplicado ao transdutor em volts (voltage), atraso de tempo do sinalcoletado em microsegundos (delay) e taxa de amostragem do sinal (rate).

Tabela 4.1: Parametros de aquisicao de dados do programa SADUS.

Parametro Descricao VariacaoBuffer tamanho do vetor em bytes 256/512/1k/16kGain ganho aplicado ao sinal recebido 1/2/5/10/20/50/100Voltage tensao do pulso no transdutor 55/90/125/195/230/265/300 VDelay atraso aplicado ao sinal amostrado 0 a 100 µsRate frequencia de amostragem 50 ou 100 MHz

A interface do programa e apresentada na Figura 4.8, [42]. A partir da versao 5.0,o programa contava com um sistema de aquisicao e apresentacao do sinal de ultra-som

64


em tempo real. Essa modificacao foi realizada em um projeto de iniciacao cientıfica degraduacao do Departamento de Engenharia Eletrica.

Figura 4.8: Interface do programa de controle e aquisicao de dados, SADUS, atualmentena versao 5.0.

Alguns outros parametros das aquisicoes tambem podem ser ajustados nas opcoesapresentadas na janela Save, como mostra a Figura 4.8. Modifica-se o numero de sinaisutilizados para se obter um sinal medio por posicao do transdutor na opcao Aqq Mean,essa pratica elimina um pouco do ruıdo aleatorio do sinal, podendo-se variar de 1a 10 aquisicoes por ponto. Tambem, o numero de posicoes para a varredura linearda amostra, Aqq, variando de 1 a 100 posicoes; o numero de passos do motor entrecada aquisicao, Steps by Aqq, variando de 1 a 10, sendo que cada passo corresponde a50 µm aproximadamente, e, por ultimo, pode-se escolher a direcao da varredura, MoveDirection.

A janela Move possui alguns controles da movimentacao do transdutor no eixolinear, utilizados somente para ajustes sem realizar medidas. Mostram-se a direcao domovimento (Direction), um botao Home que posiciona o transdutor no ponto inicialdo eixo (utilizando um sensor de fim de curso para indicacao) e dois botoes de pontosde referencia, Ref1 e Ref2, ainda nao implementados no hardware, mais o botao desaıda Quit.

Os sinais de eco coletados pelo hardware sao discretizados em amplitude usando-se 8 bits, o comprimento do sinal e definido pelo software, atraves do parametro buffer,que possui quatro possıveis valores: 256, 512, 1k e 16k, que equivalem a discretizacaoem 8, 9, 10 e 14 bits. Para se obter melhores resolucoes, utilizava-se a ultima opcao, aqual possibilita 16384 pontos para o sinal. O software salva os dados em um arquivoformato MAT, que e compatıvel com o programa Matlab, da empresa MathWorks, Inc.,

65


sendo acessıveis por protocolos em linguagem C.Quando e realizada uma varredura na amostra, repete-se esse procedimento para

cada posicao, salvando os resultados em um unico arquivo como uma matriz, sendo onumero de linhas igual ao tamanho do buffer de dados e o de colunas igual a numerode posicoes. A Figura 4.9 mostra o modelo geral dos arquivos de dados. Pode-se ver,tambem, o cabecalho usado para identificacao dos dados acrescentado ao arquivo comoa primeira coluna da matriz [42].

Figura 4.9: Modelo geral utilizado pelo software para salvar os dados em arquivos.Apresenta-se tambem o cabecalho usado para identificacao dos dados [42].

Figura 4.10: Codificacao usada no cabecalho dos arquivos para identificacao dos dados.A tensao e dada em volts, a frequencia em hertz.

66

4.3 Softwares de Simulacao e Analise

Utiliza-se somente as quatro primeiras posicoes dessa coluna, as quais conteminformacoes sobre os valores usados de ganho, tensao do pulso, taxa de amostrageme atraso de tempo inicial, conforme codificacao apresentada na Figura 4.10. O atrasoinicial nao aparece na figura, pois ele nao e codificado, colocando-se o valor em micro-segundos diretamente no cabecalho.


Os softwares de simulacao e analise foram escritos em linguagem C ANSI utilizandoa biblioteca cientıfica GSL (GNU Scientific Library), a qual e uma colecao de rotinaspara computacao numerica [10]. A programacao foi feita em ambiente RedHat Linux9.0.

O algoritmo era dividido em duas partes: simulacao e analise, sendo estas orga-nizadas da seguinte forma:

I) Simulacao

1. Forneciam-se os parametros: Ad (amplitude da parte difusa), V (jitter) e NS (numerode simulacoes).

2. Com esses valores, simulava-se o sinal de eco com o MSS desejado.

3. Forneciam-se os parametros: p (ordem da decomposicao AR) e N (tamanho da janela).

4. Janelava-se o eco e encontrava-se a envoltoria com a transformada de Hilbert.

5. Utilizava-se o metodo LMS ou da Autocorrelacao para se encontrar os coeficientes daregressao.

6. Achavam-se os polos da regressao.

7. Organizavam-se os polos pela fase.

8. Salvavam-se os resultados em arquivos.

II) Analise

1. Abriam-se os arquivos em ordem.

2. Buscavam-se, para cada sinal, os tres polos de fases positivas com maior modulo.

3. Comparava-se com o valor de MSS esperado para cada uma das harmonicas e calculava-se o erro de estimativa.

No restante desta secao, serao apresentadas as implementacoes de cada uma des-sas etapas. Na secao seguinte, alguns aspectos computacionais e metodos numericosmais especıficos serao melhor discutidos.

A parte da Simulacao possui oito etapas, sendo a que consome o maior tempocomputacional do processo. Na primeira, passavam-se os parametros que seriam usadospara a simulacao dos ecos: amplitude percentual da parte difusa do eco em relacao aparte regular (Ad), jitter percentual da parte regular (V ) e numero de simulacoes (NS)

67


para cada conjunto de valores de Ad e V . Os valores utilizados para Ad foram entre 1e 55% e, para V , entre 1 e 45%. O valor tıpico para N foi de 512 simulacoes.

A segunda etapa e de extrema importancia pois e nesta que se simulam os sinaisde eco para estudo. O gerador de ecos usado foi retirado do trabalho [22], sendo umaaplicacao dos modelos apresentados no Capıtulo 2.

Esse gerador simula os ecos retroespalhados pelo meio pela convolucao da repostaimpulsiva do modelo do tecido com o pulso de ultra-som. O modelo utilizado apre-senta quantidade de partıculas difusas igual a 7, 5 partıculas por milımetro, sendo avelocidade do som igual a 1540 m/s (tecido mole) e taxa de amostragem de 25 MHz.Os parametros do pulso de ultra-som foram: frequencia central de 2, 5 MHz, largurade banda de 2 MHz e perıodo de repeticao do pulso de 1 µs [22].

O sinal de eco e simulado utilizando-se 1000 partıculas difusas em 135 mm e aspartıculas regulares em 100 mm. O modelo do tecido e baseado na distribuicao daspartıculas periodicas por uma funcao de densidade de probabilidade gamma e as difusaspor uma funcao uniforme. A relacao sinal-ruıdo SNR do eco e calculada conformeapresentado em [22]. Para ilustracao do algoritmo gerador de ecos, apresenta-se umexemplo, considerando-se MSS de 1, 25 mm, Ad = 15% e V = 10%, como mostraa Figura 4.11, sendo (a) o pulso de ultra-som, (b) o PSD do pulso, (c) a respostaimpulsiva do meio e (d) o eco simulado.

Figura 4.11: Ilustracao dos resultados do algoritmo gerador de ecos. (a) Pulso deultra-som. (b) PSD do pulso. (c) Resposta impulsiva do meio. (d) Eco simulado.

68


A terceira etapa consiste em fornecer os parametros: ordem da regressao AR (p)e tamanho da janela (N) a ser usada no eco para analise. Segundo [44], a ordem daregressao usada para estimar o MSS nao pode exceder N/5 e deve ser em torno de25% maior (ou um pouco mais) que o numero de amostras usadas pelo MSS, ou seja,MMSfs/c, sendo fs a frequencia de amostragem e c a velocidade do som no modelo(1540 m/s). Dessa forma, como o tamanho da janela usado foi 512 pontos, e os valoresde MSS variavam entre 1, 25 e 3 mm, utilizou-se uma ordem igual a 64.

Na quarta etapa, encontra-se o sinal janelado e obtem-se sua envoltoria pelo valorabsoluto de sua transformada de Hilbert, conforme sera apresentado na proxima secao.Em seguida, retira-se o valor medio do sinal.

A quinta etapa consiste em encontrar os coeficientes da decomposicao AR deordem p do sinal do item anterior. Para isso, dois metodos foram utilizados, primei-ramente o algoritmo LMS e, por fim, o metodo da autocorrelacao. O uso do LMSbaseia-se na equivalencia do modelamento all-pole e da predicao linear como apresen-tado no Capıtulo 3. Buscavam-se os coeficientes que minimizassem o erro entre y[n]e sua estimativa y[n], composta por uma combinacao linear de p valores passados,y[n− k], com k = 1, . . . , p.

Entretanto, essa abordagem, descrita com mais detalhes na secao seguinte, apre-sentava erros nos polos encontrados na sexta etapa, portanto foi substituıda pelometodo da autocorrelacao. Os detalhes de programacao serao discutidos na proximasecao.

Na sexta etapa, encontravam-se os polos da decomposicao AR utilizando-se ummetodo de calculo de raızes de polinomios, no caso o polinomio de grau p formado pelosak encontrados, sendo k = 0, . . . , p. Os detalhes de programacao serao apresentados naproxima secao. Na setima e oitava etapas, organizavam-se os polos pelos seus valoresde fase e salvavam-se os resultados em arquivos para serem usados na parte seguinte.

A parte da analise e composta por tres etapas, com pouco custo computacio-nal. A primeira consiste somente em abrir sequencialmente os arquivos salvos comas simulacoes. Para cada combinacao de Ad e V , escolhia-se um numero que seriausado para nomear os arquivos, entao, na analise, abriam-se os arquivos na ordem danumeracao.

Na segunda etapa, conforme apresentado no Capıtulo 3, tinha-se o objetivo deencontrar os tres polos com maior modulo, os quais corresponderiam aos picos de maioramplitude das tres primeiras harmonicas do espectro de potencia do sinal. Para isso,analisavam-se somente os polos com fases positivas, pois os outros eram conjugados.Entao, achava-se o maximo dos modulos, zerando o seu valor em seguida, para repetiro processo e encontrar o segundo e terceiro maiores polos. Por ultimo, calculava-se oerro de estivativa dos valores de MSS encontrados com cada um dos polos.

Para ilustrar os procedimentos apresentados nesta secao, utiliza-se um exemplo deestimativa das frequencias de uma soma de duas senoides com o algoritmo proposto. Asfrequencias sao 50 e 75 Hz, entao, segundo [2], para rastrear um sinal composto de duassenoides, deve-se utilizar um modelo autorregressivo de quarta ordem (o dobro), p = 4.Foi usada uma frequencia de amostragem de 162Hz, que respeita a taxa de Nyquist.Os polos da decomposicao encontrados foram: −0, 8583 + 0, 2064i; −0, 3416 + 0, 9296ie seus conjugados, assim, as frequencias das senoides sao 49, 83 e 75, 29 Hz, um erro de0, 34 e 0, 39% respectivamente. Neste exemplo, nao se usa a envoltoria, pois o objetivo

69


e rastrear os senos.

Para uma melhor ilustracao, apresenta-se um exemplo usando o sinal de eco daFigura 4.11 (d), sendo MSS igual a 1, 25 mm, Ad = 15% e V = 10%. Para esse casoo valor da frequencia esperada e fMSS = 0, 616 MHz. A Figura 4.12 apresenta aenvoltoria da janela do eco em (a) e o PSD desse sinal em (b), onde se pode ver ospicos no espectro. A Figura 4.13 apresenta os modulos dos polos da decomposicao ARem funcao da frequencia (ate 2 MHz) utilizando os dois metodos, (a) autocorrelacaoe (b) LMS. Os erros dos tres primeiros harmonicos para esses dois casos sao: −2, 74;−2, 82 e −4, 39% para o metodo da autocorrelacao e −8, 85; −2, 45 e 1, 52% para oLMS. Pode-se ver que, para o LMS, o erro para a primeira harmonica e quatro vezesmaior que o erro da autocorrelacao.

Figura 4.12: (a) Envoltoria de uma janela de 512 pontos do eco apresentado na Figura4.11, com Ad = 15% e V = 10%. (b) PSD desse sinal janelado mostrando os picos dastres primeiras harmonicas.

4.4 Aspectos Computacionais

Na secao anterior, apresentaram-se os procedimentos usados para simulacao e analisede sinais de ultra-som, mas alguns aspectos computacionais foram omitidos para maiorobjetividade da explanacao. Assim, nesta secao, esses pontos sao retomados com maisprofundidade, adentrando os metodos numericos utilizados, os quais foram baseadosprincipalmente no livro Numerical Recipes in C [34].

Dos algoritmos apresentados anteriormente, quatro precisaram de metodos refi-nados para serem desenvolvidos, estes sao: calculo da convolucao de dois sinais, calculoda transformada de Hilbert de um sinal, os metodos da autocorrelacao e LMS, e, por

70


Figura 4.13: Graficos com os modulos dos polos da decomposicao AR em funcao dafrequencia (ate 2 MHz) utilizando os dois metodos propostos. (a) Metodo da autocor-relacao. (b) Algoritmo LMS.

ultimo e mais complicado, calculo das raızes complexas de um polinomio de alta ordeme de autovalores de matrizes. Esses quatro pontos serao explicados, um a um, a seguir.

Convolucao

O calculo da convolucao entre dois sinais era importante pois, os ecos retroes-palhados simulados sao obtidos pela convolucao do modelo impulsivo do meio com opulso de ultra-som.

Foram utilizadas duas abordagens para esse calculo, uma direta e mais lenta eoutra indireta e mais rapida. A primeira, utiliza a definicao de convolucao [30], comomostra a Equacao 4.3, sendo y[n], o resultado da convolucao de h[n] e x[n].

y[n] =∞∑

k=−∞

x[k]h[n− k] (4.3)

Mas, na pratica, os sinais x e h nao sao infinitos, mas sao definidos em intervalos[0, Nx] e [0, Nh] respectivamente. Entao, um artifıcio utilizado para resolver esse pro-blema foi o uso de matrizes. Assim, a equacao anterior transforma-se na Equacao 4.4,sendo ~x[n] = x[0], . . . , x[Nx − 1], ~H, de dimensoes Nx por Nx + Nh + 1, apresentadoa seguir e o vetor convoluıdo ~y[n] = y[0], . . . , y[Nx +Nh], com dimensao Nx +Nh +1.

~y[n] = ~x[n] · ~H[n] (4.4)

~H[n] =

h[0] . . . h[Nh − 1] 0 0 . . . 00 h[0] . . . h[Nh − 1] 0 . . . 0...

.... . .

......

. . ....

0 0 . . . 0 h[0] . . . h[Nh − 1]

71


Para melhor compreensao, apresenta-se um exemplo bem simples, onde ~x = [1, 2]

e ~h = [3, 4, 5], assim, tem-se ~y[n] = y[0], y[1], y[2], y[3]. Expandindo-se a Equacao 4.3tem-se as relacoes:

y[0] = x[0]h[0]+ x[1]h[−1]y[1] = x[0]h[1]+ x[1]h[0]y[2] = x[0]h[2]+ x[1]h[1]y[3] = x[0]h[3]+ x[1]h[2]

Assim, usando-se a Equacao 4.4, encontra-se as relacoes anteriores no formatomatricial apresentado a seguir, sendo o resultado ~y[n] = [3, 10, 13, 10].

y[0] y[1] y[2] y[3]

=

x[0] x[1]·

h[0] h[1] h[2] 00 h[0] h[1] h[2]

Entretanto, essa abordagem direta para calculo da convolucao e computacional-

mente custosa, pois, para sinais com muitos valores, a multiplicacao matricial consomemuito tempo. Dessa forma, buscou-se um outro metodo, este de forma indireta, porqueutiliza a transformada de Fourier. Ele consiste em achar as transformadas dos sinais,multiplicar esses valores e achar a transformada inversa, conforme apresentado pelaEquacao 4.5, sendo F o operador transformada de Fourier.

Y (ejw) = Fh[n] ∗ x[n] = H(ejw)X(ejw) ⇒ y[n] = F−1H(ejw)X(ejw) (4.5)

Com esse metodo indireto, pode-se calcular a convolucao utilizandos os rapidosalgoritmos de FFT (Fast Fourier Transform) Radix-2 que a biblioteca GSL [10] pos-sui. Esse algoritmo considera que os sinais tenham dimensao potencia de dois, entaopreenchem-se os vetores com zero (zero padding) ate a dimensao potencia de dois maisproxima.

Transformada de Hilbert

A transformada de Hilbert e uma transformacao que relaciona as partes reaise imaginarias de uma transformada de Fourier (TF), essa relacao e unica quando asequencia em questao e causal [30]. O conceito por tras dessa transformada e quequalquer sequencia x[n] pode ser escrita como a soma de uma sequencia par, xe[n],e outra ımpar, xo[n], podendo-se expressar esses termos conforme apresentado pelaEquacao 4.6.

xe[n] =x[n] + x[−n]

2e xo[n] =

x[n]− x[−n]

2(4.6)

Assim, caso a sequencia x[n] seja causal, ou seja, x[n] = 0 para n < 0, pode-serecuperar essa sequencia a partir da parte par ou da parte ımpar, sendo que, nesteultimo caso, somente para n 6= 0, conforme as Equacoes 4.7 e 4.8.

x[n] = 2xe[n]u[n]− xe[0]δ[n] (4.7)

x[n] = 2xo[n]u[n] + x[0]δ[n] (4.8)

72


Entao, caso x[n] seja estavel, a sua transformada de Fourier, X(ejw), existe e,se x[n] for real, a TF de xe[n] e a parte real de X(ejw), XR(ejw), e a TF de xo[n] ea parte imaginaria XI(e

jw). Dessa forma, a parte real da transformada da sequenciaa determina completamente, e XI(e

jw) tambem pode ser determinada pela parte real.Conforme [30], essas relacoes, chamadas de transformacoes de Hilbert, podem ser esci-tas conforme as Equacoes 4.9 e 4.10, sendo P o valor principal de Cauchy.

XI(ejw) = − 1

2πP∫ π

−π

XR(ejθ)cot(w − θ

2)dθ (4.9)

XR(ejw) = x[0] +1

2πP∫ π

−π

XI(ejθ)cot(

w − θ

2)dθ (4.10)

Percebe-se que um algoritmo de calculo da transformada de Hilbert (TH) poressas equacoes seria algo muito complicado, assim, um artifıcio deve ser usado. Nocaso, o que se deseja fazer e encontrar X(ejw) a partir de XR(ejw), sendo necessarioutilizar a Equacao 4.9, mas, segundo [27], a TH de uma funcao real e estavel x(t),definida em −∞ < t < ∞, e equivalente a convolucao entre essa funcao e 1/πt, assim,a TF da parte imaginaria de uma funcao analıtica z(t) que corresponde a x(t) e definidapela Equacao 4.11.

ZI(f) =

X(f), para f > 0

0, para f = 0−X∗(−f), para f < 0

(4.11)

Assim, considerando o caso discreto, x[n], uma funcao real de tamanho (N) par(caso seja ımpar, basta somente usar N = N + 1), tem-se que z[n], uma funcao tipoanalıtica [27] correspondendo a x[n], pode ser calculada utilizando-se as Equacoes 4.12e 4.13, sendo Z[m] e X[m] as TF de z[n] e x[n] respectivamente e Ts o perıodo deamostragem. Assim, mais uma vez, utiliza-se o algoritmo FFT Radix-2 para calcularas transformadas de Fourier e suas inversas.

Z[m] =

X[0], para m = 0

2X[m], para 1 ≤ m ≤ N/2− 1X[N/2], para m = N/2

0, para N/2 + 1 ≤ m ≤ N − 1

(4.12)

zI(f) =1

NTs

N−1∑m=0

Z[m]ej2πmn/N (4.13)

Metodos da Autocorrelacao e LMS

Os algoritmos que executam os metodos da Autocorrelacao e LMS foram pro-gramados usando as formas apresentadas em [13]. Ambos dependem de uma matriz

de convolucao, ~X, obtida do sinal de entrada, x[n], de tamanho N . Considerando aordem escolhida p, essa a matriz possui dimensoes N + p− 1 por p.

73


As colunas de ~X consistem em deslocamentos do sinal de entrada, sendo que nacoluna i, o sinal e x[n− i], preenchendo com zeros os termos desconhecidos. A Equacao4.14 apresenta uma forma geral da matriz de convolucao do sinal x[n].

~X =

x[0] 0 . . . 0 0... x[0] . . .

......

x[N − 1]...

. . .... 0

0 x[N − 1] . . . x[0] 0

0 0. . .

... x[0]...

... . . . x[N − 1]...

0 0 . . . 0 x[N − 1]

(4.14)

Assim, utilizando-se as equacoes dos filtros LMS apresentadas no Capıtulo 3, 3.54e 3.60 e lembrando-se que e[n] = d[n]− d[n], obtem-se o procedimento a seguir, sendolini., o operador que seleciona a linha i de uma matriz:

1. e[0] = d[0]

2. ~w[0] = µ e[0] lin0 ~X∗

3. para 1 ≤ k ≤ N − 1

e[k] = d[k]− ~w[k − 1] link ~XT

~w[k] = ~w[k − 1] + µ e[k] link ~X∗

O algoritmo do metodo da autocorrelacao, ACM (Autocorrelation Method), eimplementado utilizando a Equacao 3.81 e o erro, E , e calculado pelo valor absoluto doproduto interno de ~x1 com −~x1, como pode ser visto no procedimento a seguir, sendored. um operador de reducao de ordem, que aplicado a uma matriz ~A de dimensoes Mpor N , retorna uma matriz com as linhas 0 a M−2 e colunas 0 a N−2 da matriz original.Tem-se que ~a1 = a[1], a[2], . . . , a[p], ~x1 = X[1, 0], X[2, 0], . . . , X[N + p − 1, 0],pinv. e o operador pseudoinversa, a ser explicado mais a frente, e abs. e o operadorvalor absoluto.

1. ~Xp−1 = red ~X

2. a[0] = 1

3. ~a1 = −pinv ~Xp−1 ~x1

4. E = abslin0 ~XT ~X ~a

O algoritmo apresentado necessita do calculo da pseudoinversa de uma matriz,o que nao e simples. A abordagem utilizada baseia-se na aplicacao da decomposicaoem valores singulares, SVD (Singular Value Decomposition), da matriz ~Xp−1. Essadecomposicao, implementada na biblioteca GSL [10], pode ser aplicada a qualquer

matriz ~A, M por N , resultando em um produto da matriz ortogonal ~U , M por N , da

74


matriz diagonal ~S de valores singulares, N por N , e a transposta da matriz ortogonal~V , N por N , [12].

Uma matriz ortogonal e aquela cuja inversa e a sua propria transposta, ou seja,o produto dela por sua transposta e igual a identidade. A Equacao 4.15 apresenta adecomposicao da matriz ~A, sendo os valores singulares σi = Sii nao-negativos. Apos aSVD, pode-se encontrar a pseudoinversa da matriz ~A utilizando a Equacao 4.16, se orank for igual a N , ou seja, nenhum valor singular e zero, se esse for o caso, deve-seentao, usar somente os R valores nao-nulos, transformando ~S em uma matriz quadrada,R por R, e ~U e ~V da mesma forma.

~A = ~U ~S ~V T (4.15)

pinv ~A = ~V (~ST ~S)−1 ~ST ~UT (4.16)

Raızes de Polinomios e Autovalores

O algoritmo de calculo das raızes complexas de um polinomio de alta ordemfoi a parte mais complicada da implementacao dos softwares. Tentou-se utilizar osmetodos da linearizacao e de Bairstow (Capıtulo 3), mas estes apresentaram problemasde convergencia e de escolha de valores iniciais. Dessa forma, fez-se necessario escolherum outro metodo, o qual, mesmo que mais lento, fornecesse resultados confiaveis.

Utiliza-se a solucao pela matriz companheira [33, 34], pois, sabe-se que as raızesde um polinomio de grau p, P (x) = a0x

p + a1xp−1 + . . . + ap−1x + ap, sao iguais aos

autovalores da matriz companheira, ~C, com dimensao p por p, definida pela Equacao4.17.

~C =

−a1

a0

−a2

a0

. . . −ap−1

a0

−ap

a0

1 0 . . . 0 00 1 . . . 0 0...

.... . .

......

0 0 . . . 1 0

(4.17)

O problema e que as bibliotecas de calculo de autovalores do GSL so se aplicama matrizes reais simetricas (matrizes que sao iguais a suas transpostas) ou complexashermitianas. Esse nao e o caso da matriz companheira, que e Hessenberg superior, poispossui os elementos c(i, j) = 0 para i > j + 1, ou seja, so tem zeros abaixo da primeirasubdiagonal inferior [12].

Entao, para se resolver esse problema, implementou-se um algoritmo de calculode autovalores (sem os autovetores) para matrizes Hessenberg utilizando decomposicaoQR com deslocamento simples para evitar operacoes complexas [43]. O algoritmo dedecomposicao QR usado faz parte das bibliotecas GSL.

O conceito por tras da decomposicao QR e que qualquer matriz real, ~A, pode serdecomposta em um produto de duas matrizes ~Q e ~R, sendo a primeira, ortogonal, e asegunda, triangular superior. Considerando a matriz, ~A′, formada pelo produto ~R ~Q,e lembrando que ~Q−1 = ~QT devido a propriedade ortogonal, tem-se a Equacao 4.18.

75


~A′ = ~QT ~A ~Q (4.18)

Assim, ~A′ e uma transformacao ortogonal de ~A, entao, tambem possui o formatoHessenberg [34] e os mesmos autovalores. A forma de se encontrar esses autovalores

se baseia em repetidas transformacoes QR ate que ~A se transforme em triangularinferior com os autovalores nas diagonais, caso sejam todos reais e diferentes. Nocaso de multiplicidade m, a matriz converge para triangular, exceto em um blocomatricial diagonal de ordem m, cujos autovalores sao os de multipla ordem. Caso osautovalores sejam complexos, ocorre algo parecido, criam-se blocos de ordem dois, comos autovalores. Para diminuir o numero de calculos, utiliza-se a deflacao da matriz,[34], ou seja, ir reduzindo a sua ordem conforme se encontram os autovalores.

Uma tecnica usada para acelerar o processo e o deslocamento, pois se s e umaconstante, entao os autovalores da matriz ~A−s~I sao λi−s, sendo ~I a matriz identidadee λi os autovalores de ~A. Dessa forma, tem-se o processo a seguir, sendo o ındice k onumero da iteracao da recursao:

~Ak − sk~I = ~Qk

~Rk

~Ak+1 = ~Rk~Qk + sk

~I

O deslocamento e escolhido como o ultimo elemento da diagonal principal damatriz (canto inferior direito). A iteracao acaba quando c[p − 1, p − 2] = 0, mas, napratica, utilizou-se c[p − 1, p − 2] < 0, 5.10−4, dessa forma, tem-se que c[p − 1, p − 1]e um autovalor (real). Em seguida, retiram-se a linha e coluna p − 1 da matriz e serecomeca o procedimento. Outra condicao de parada e quando c[p − 2, p − 3] = 0,entao os autovalores da submatriz 2 por 2 do canto inferior direito sao os autovaloresprocurados, que serao complexos conjugados. Assim, apagam-se a linha e coluna p− 1e p− 2, recomecando o processo.

Apresentam-se, a seguir, exemplos para ilustracao desse procedimento. Considera-se a matriz ~A a seguir, a qual, apos quatro iteracoes, apresenta um zero (valor menorque 0, 5.10−4) na posicao (4, 3), como se pode ver, sendo−0, 638 um de seus autovalores.

~A =

4 0, 5 1 −1, 5 −21 0 0 0 00 1 0 0 00 0 1 0 00 0 0 1 0

⇒

4, 153 −1, 899 0, 936 −1, 515 0, 2480, 004 0, 461 0, 176 0, 875 −0, 338

0 0, 602 0, 139 −0, 580 −0, 4110 0 0, 942 −0, 115 0, 3080 0 0 0 −0, 638

No exemplo a seguir, encontra-se um par de autovalores complexos. Apos 13

iteracoes, chega-se a uma matriz onde a[2, 1] = 0, entao, a partir da matriz 2 por2 formada no canto inferior direito, encontram-se os autovalores, −0, 449 + 0, 480i e−0, 449− 0, 480i.

~A =

4 0, 5 −1 −1, 51 0 0 00 1 0 00 0 1 0

⇒

4, 040 −1, 951 −0, 286 0, 2160 0, 859 0, 098 0, 5540 0 −0, 679 0, 3890 0 −0, 728 −0, 219

76


[−0, 679 0, 389−0, 728 −0, 219

]⇒ (λ + 0, 679)(λ + 0, 219)− (−0, 728 × 0, 389) = 0

Mas, como se pode ver nesse exemplo, o elemento a[1, 0] tambem e zero, sendoque, na verdade, isso aconteceu em torno da nona iteracao. Isso e algo muito comumde se acontecer, principalmente com matrizes grandes, entao, utiliza-se outros testesem conjunto com os ja citados para melhorar a performance do algoritmo. Testa-sese algum outro elemento da subdiagonal inferior e zero, se isso acontecer, armazena-sea parte da matriz acima desse ponto e se trabalha com a parte de baixo. Apos seresolver essa parte, volta-se para a que tinha sido armazenada. Esse processo pode serfeito varias vezes, assim, o custo computacional do algoritmo diminui. O procedimentoa seguir resume o algoritmo de calculo de raızes de polinomios utilizado, sendo δ =0, 5.10−4:

1. Encontra-se a matriz companheira, ~C, a partir dos coeficientes;

2. enquanto p > 0 faca:

(i) se p = 1 faca:

- o autovalor e c[1, 1]

- faz-se p = 0, encerrando a etapa

(ii) se p = 2 faca:

- os autovalores sao os autovalores da submatriz inferior direita

- faz-se p = 0, encerrando a etapa

(iii) se p > 2 faca:

(a) sk = c[p− 1, p− 1]

(b) aplica-se a decomposicao QR para ~Ck − sk~I

(c) encontra-se ~Ck+1 = ~Rk~Qk + sk

~I

(iv) se c[p− 1, p− 2] < δ, um autovalor e c[p− 1, p− 1]

- decrementa p em uma unidade

- retira-se a linha e coluna p− 1

(v) se c[p−2, p−3] < δ, dois autovalores sao os autovalores da submatriz inferiordireita

- decrementa p em duas unidades

- retira-se as linhas e colunas p− 1 e p− 2

(vi) se c[m, m− 1] < δ com m < p− 2, divide-se a matriz, armazena-se a partesuperior e continua na inferior

- decrementa p em m− 1 unidades

3. se alguma matriz foi armazenada, recupera-la e voltar para a etapa anterior

Para ilustrar esse procedimento, apresenta-se um exemplo: encontrar os polos dafuncao de sistema de quinta ordem: H(z) = 1/(1−1, 6z−1−0, 2z−2 +0.4z−3 +0, 6z−4 +

77


z−5). Primeiramente, escreve-se a matriz companheira, ~C0, lembrando que o polinomioe P (z) = z5 − 1, 6z4 − 0, 2z3 + 0.4z2 + 0, 6z + 1, assim, tem-se:

~C0 =

1, 6 0, 2 −0, 4 −0, 6 −11 0 0 0 00 1 0 0 00 0 1 0 00 0 0 1 0

Em seguida, sendo p = 5, escolhe-se s0 = c[4, 4], aplica-se a decomposicao QR e

encontra-se ~C1 = ~R0~Q0 + s0

~I, como apresentado a seguir:

~C1 =

1, 690 0, 247 −0, 681 1, 036 −0, 478−0, 533 −0, 112 −0, 142 0, 193 −0, 814

0 −1, 016 0, 086 −0, 131 0, 0560 0 −1, 024 0, 082 0, 2510 0 0 −0, 471 −0, 146

Como se pode ver, nenhum dos elementos da primeira subdiagonal inferior,

−0, 533; −1, 016; −1, 024 e −0, 471, sao menores que δ = 0, 5.10−4. Entao, continua-sea recursao, buscando-se ~C2 pela decomposicao QR. Para simplicacao, considera-se aquinta iteracao, k = 5, quando a matriz ~C5 e:

~C5 =

1, 383 1, 408 −0, 146 −0, 114 0, 0724−0, 138 1, 369 0, 193 0, 0236 −0, 116

0 0, 127 −0, 331 0, 908 0, 1510 0 −0, 663 −0, 0433 −0, 4700 0 0 0, 0001 −0, 778

Entao, tem-se c[4, 3] < δ, assim, λ0 = c[4, 4] = −0, 778. Para continuar, decrementa-

se a ordem, p = 4 e retira-se a linha e coluna 4, reduzindo-se a matriz. Volta-se para asegunda etapa do procedimento, buscando ~C6. Na decima primeira iteracao (k = 11),

encontra-se a seguinte matriz, ~C11:

~C11 =

1, 561 0, 164 −0, 090 −0, 084−1, 374 1, 201 −0, 256 −0, 044

0 −0, 0004 −0, 291 0, 9500 0 −0, 615 −0, 094

Assim, tem-se c[2, 1] < δ, entao, encontrou-se dois autovalores complexos conju-

gados, os valores sao obtidos da submatriz 2 por 2 do canto inferior direito:∣∣∣∣ λ + 0, 291 −0, 950+0, 615 λ + 0, 094

∣∣∣∣ = 0 ⇒ λ2 + 0, 385λ + 0, 612 = 0

Dessa forma, os autovalores sao λ1 = −0, 192 + 0, 758i e λ2 = −0, 192 − 0, 758i.Decrementa-se a ordem em duas unidades, p = 2, e retira-se as linhas e colunas 3 e4. Como a matriz resultante e 2 por 2, nao e necessario fazer a decomposicao QR,bastando encontrar os autovalores da matriz que sobrou:

78


∣∣∣∣ λ− 1, 561 −0, 164+1, 374 λ− 1, 201

∣∣∣∣ = 0 ⇒ λ2 − 2, 763λ + 2, 101 = 0

Assim, os dois ultimos autovalores sao λ3 = 1, 381+0, 439i e λ4 = 1, 381−0, 439i.Dessa forma, os polos de H(z) sao esses autovalores encontrados.

79

80

Capıtulo 5

Resultados

Neste capıtulo serao apresentados os resultados obtidos pela aplicacao da decomposicaoautorregressiva e rastreamento de polos ao sinais retroespalhados de ultra-som criadosa partir dos modelos dos meios. Assim, pode-se interpretar melhor os fenomenos en-volvidos nessas analises.

Em seguida, baseando-se nessa compreensao das acoes do metodo sobre os si-nais, serao apresentados os resultados de simulacoes de Monte Carlo, variando-se osparametros de regularidade do meio (jitter) e relacao sinal-ruıdo (SNR). Assim, pode-se avaliar os efeitos que essas variacoes, as quais estao muitas vezes relacionadas comdisfuncoes ou doencas nos tecidos, acarretam nos resultados.

Por ultimo, serao apresentados os resultados da aplicacao das tecnicas em sinaisretroespalhados reais, coletados de phantoms feitos com gelatina e po de grafite, osquais simulam tecidos biologicos moles. Dessa forma, uma compreensao ainda maiorpode ser conseguida, pois esses sinais possuem efeitos de caracterısticas nao utilizadasnos modelos propostos.

5.1 Resultados com Sinais Artificiais

Nesta secao, apresentam-se os resultados para analises feitas aplicando os metodos asinais artificiais. Para uma boa compreensao dos efeitos da analise, utilizam-se quatrocasos de teste. A Tabela 5.1 mostra alguns parametros que foram mantidos fixos nostestes. Como visto, o MSS foi mantido constante, assim, utilizando a Equacao 2.7encontra-se o valor de frequencia do pico do espectro referente a esse espacamentoigual a fMSS = 0, 308 MHz.

Os parametros que foram variados para as analises sao a regularidade (jitter)e a relacao sinal-ruıdo (SNR), sendo este um parametro relacionado aos processosaleatorios usados para gerar o eco, mas dependente da amplitude percentual do ruıdoutilizada, Ad, como apresentado na Tabela 5.2.

Utilizando-se os algoritmos descritos no Capıtulo 4, analisaram-se sinais de ecoartificiais construıdos com esses parametros. Nas Figuras 5.1 a 5.4, apresentam-se,respectivamente para os casos de numero 1 a 4, os sinais de envoltoria com media nula(em cima) junto com seus espectros de potencia (embaixo).

81

Capıtulo 5 - Resultados

Tabela 5.1: Parametros mantidos constantes nas analises dos quatro casos de teste.

Parametro Valor

Velocidade do som no meio 1540 m/sFrequencia de amostragem 25 MHzFrequencia de ressonancia do transdutor 2, 5 MHzLargura de banda do transdutor 2, 5 MHzTamanho da janela 512 amostrasOrdem da decomposicao AR 64Espacamento medio dos espalhadores 2, 5 mm

Tabela 5.2: Parametros variados nas analises dos quatro casos de teste.

Numero Caso de Teste Jitter SNR Ad

Caso 1 Alta regularidade e baixo ruıdo 5% 7,894 dB 20%Caso 2 Alta regularidade e alto ruıdo 5% 1,120 dB 40%Caso 3 Baixa regularidade e baixo ruıdo 25% 6,904 dB 20%Caso 4 Baixa regularidade e alto ruıdo 25% 1,321 dB 40%

Figura 5.1: Sinal de envoltoria com media nula usado na analise do caso de teste 1 (emcima) e espectro de potencia do sinal (embaixo).

82




83



O metodo usado consiste no rastreamento dos polos da decomposicao AR do sinal,lembrando que cada par de polos complexos gera um pico no espectro de frequencias, oqual possui amplitude maior quanto mais perto do cırculo unitario estiverem os poloscorrespondentes.

Para se encontrar a frequencia de MSS, fMSS, busca-se o polo positivo com maiormodulo, pois este esta associado ao primeiro pico no espectro de potencia. Entao,para melhor visualizacao, as Figuras 5.5 a 5.8 apresentam os modulos dos polos dadecomposicao em funcao de sua frequencia, respectivamente para os casos de numero 1a 4. Nessas figuras, somente aparecem os polos com frequencias positivas, pois a partenegativa e simetrica.

Os resultados das analises para os quatro casos de teste sao apresentados naTabela 5.3. Percebe-se que, para alta regularidade e baixo ruıdo, obteve-se o menorerro, para os casos 2 e 3, os erros foram intermediarios e, para o ultimo caso, baixaregularidade e alto ruıdo, tem-se o maior erro.

Para se encontrar valores mais fieis ao que realmente acontece com a variacaodesses parametros, deve-se utilizar uma media de varias analises, por isso, utiliza-se asimulacao de Monte Carlo na proxima secao. Entretanto, para apresentar uma visaogeral do comportamento do sistema, realizou-se uma analise baseada na media (e desviopadrao) para 100 simulacoes dos mesmos quatro casos, apresentando os resultados naTabela 5.4. A Figura 5.9 mostra uma comparacao das variacoes do erro para os quatrocasos.

84


Figura 5.5: Modulo dos polos da decomposicao AR do sinal usado no caso 1 em funcaoda frequencia (somente a parte positiva do espectro).


85




86

5.2 Simulacoes de Monte Carlo

Tabela 5.3: Resultados encontrados para os quatro casos de teste.

Numero Caso de Teste fMSS MSS Erro

Caso 1 Alta regularidade e baixo ruıdo 0,294 MHz 2,6 mm -4,63%Caso 2 Alta regularidade e alto ruıdo 0,283 MHz 2,7 mm -8,84%Caso 3 Baixa regularidade e baixo ruıdo 0,284 MHz 2,7 mm -8,28%Caso 4 Baixa regularidade e alto ruıdo 0,276 MHz 2,8 mm -11,53%

Tabela 5.4: Resultados medios encontrados em 100 simulacoes para os quatro casos deteste. Apresentam-se a media (µ) e desvio padrao (σ) dos resultados.

Caso de Teste fMSS (MHz) MSS (m/s) Erro (%)

µ σ µ σ µ σCaso 1 0,276 0,019 2,80 0,200 -12,06 8,01Caso 2 0,290 0,010 2,66 0,096 -6,47 3,83Caso 3 0,304 0,055 2,71 1,344 -8,59 53,76Caso 4 0,305 0,063 2,66 0,718 -6,32 28,72


Utilizaram-se simulacoes de Monte Carlo para se buscar uma caracterizacao do metodoquanto a sensibilidade em relacao a variacoes dos parametros: regularidade (jitter) equantidade de ruıdo (SNR). Para isto, fizeram-se quatro analises: sensibilidade aojitter, variando-se esse parametro e observando-se os resultados quanto a deteccao doMSS do sinal; sensibilidade ao ruıdo, variando-se o Ad, que e a amplitude percentualda parte difusa do modelo do eco em relacao a parte regular; sensibilidade a variacao deambos os parametros; e por ultimo, uma analise da escolha da ordem da decomposicaoAR. O valor de MSS usado foi 1,25 mm.

Analise do Efeito da Variacao da Regularidade

Para a analise da variacao do jitter, utilizaram-se quatro valores: 5%, 15%, 25%e 35% com Ad igual a 15% e 30%. Para cada uma dessas combinacoes, realizaram-se512 simulacoes. Com os resultados, construiram-se histogramas, os quais apresentamo numero de ocorrencias para cada caso em funcao da frequencia do MSS encontrado.

As Figuras 5.10 a 5.13 apresentam os resultados para a primeira harmonica, sendo(a) e (b) os histogramas para os valores de Ad iguais a 15% e 30% respectivamente. AFigura 5.14 mostra um resultado consolidado para a variacao do jitter, apresentado oserros para MSS de cada uma das tres primeiras harmonicas para cada valor de jitter eAd.

87


Figura 5.9: Variacoes dos valores de erro para os quatro casos de teste.

(a) (b)

Figura 5.10: Histograma das ocorrencias dos valores de MSS encontrados para jitterigual a 5% para a primeira harmonica. (a) Ad = 15%. (b) Ad = 30%.

88


(a) (b)


(a) (b)


89


(a) (b)


(a) (b)

Figura 5.14: Resultados consolidados para as simulacoes de variacao de jitter mos-trando os erros de MSS para as tres harmonicas. (a) Ad = 15%. (b) Ad = 30%.

90


Analise do Efeito da Variacao da Quantidade de Ruıdo

Para a analise da variacao da quantidade de ruıdo, utilizaram-se quatro valo-res: 10%, 25%, 40% e 55% com jitter igual a 10% e 30%. Para cada uma dessascombinacoes, realizaram-se 512 simulacoes, construindo-se histogramas, os quais apre-sentam o numero de ocorrencias para cada caso em funcao da frequencia do MSSencontrado.

A Figura 5.15 mostra um resultado consolidado para a variacao do ruido, apre-sentado os erros para MSS de cada uma das harmonicas para cada valor de jitter e Ad.As Figuras 5.16 a 5.19 apresentam os resultados para a primeira harmonica, sendo (a)e (b) os histogramas para os valores de Ad iguais a 15% e 30% respectivamente.

(a) (b)

Figura 5.15: Resultados consolidados para as simulacoes com variacao de quantidadede ruıdo, mostram-se os erros de MSS para as tres harmonicas. (a) Jitter = 10%. (b)Jitter = 30%.

91


(a) (b)

Figura 5.16: Histograma das ocorrencias dos valores de MSS encontrados para Ad iguala 10% para a primeira harmonica. (a) Jitter = 10%. (b) Jitter = 30%.

(a) (b)


92


(a) (b)


(a) (b)


93


Analise do Efeito de Ambas as Variacoes

Para se completar as analises de variacao de jitter e quantidade de ruıdo, realizaram-se 512 simulacoes para 25 combinacoes desses dois parametros, observando-se a quan-tidade de acertos com duas consideracoes: margem de 5% e de 25% ao redor do valoresperado, conforme [22, 35].

Os valores usados foram iguais para os dois parametros, sendo 5%, 15%, 25%,35% e 45%. Para observacao dos resultados, construiram-se figuras compostas por doisgraficos, um para margem de 5% e outro para 25%. Nestes, apresentavam-se o per-centual de acertos em funcao do jitter, utilizando-se cinco curvas, uma para cada valorde Ad. Fez-se tambem uma analise dos acertos para os primeiros tres harmonicos. AsFiguras 5.20 a 5.22 apresentam esses resultados para as harmonicas 1 a 3 respectiva-mente. Os valores de Ad usados representam uma relacao sinal-ruıdo media de 19,00;9,45; 5,05; 2,20 e -0,10 dB respectivamente.

(a) (b)

Figura 5.20: Percentual de acertos medios do MSS para a primeira harmonica em 512simulacoes variando-se o jitter e Ad. (a) Margem de 5%. (b) Margem de 25%.

Para uma visao mais geral da resposta do sistema a essas variacoes, realizaram-se outras simulacoes, mas agora, com 32 valores para cada parametro. Esses valoresforam escolhidos de forma a apresentarem crescimento exponencial. Utilizou-se a regraapresentada na Equacao 5.1, sendo p(i), um parametro (jitter ou Ad) para uma dadaiteracao i, a qual varia de 2 a 32, com valor inicial igual a 1. Os valores obtidos saoapresentados na Figura 5.23.

94


(a) (b)

Figura 5.21: Percentual de acertos medios do MSS para a segunda harmonica em 512simulacoes variando-se o jitter e Ad. (a) Margem de 5%. (b) Margem de 25%.

p(i) = p(i− 1)× 1, 125 (5.1)

Utilizando-se o erro medio dos valores de MSS obtidos para as tres harmonicasnas 512 simulacoes para cada condicao apresentada, construıram-se as Figuras 5.24 a5.26, as quais mostram graficos tridimensionais do erro em funcao dos valores de jittere Ad.

Mas, como se sabe que o valor esperado de fMSS e 0,616; 1,232 e 1,848 MHzpara cada uma das tres primeiras harmonicas respectivamente, pode-se desconsideraros polos com frequencias acima de 2 MHz. Entao, os erros muito elevados (maioresque 100%) que ocorriam devido a polos de altas frequencias com modulos proximos daunidade sao evitados. As Figuras 5.27 a 5.29 apresentam os resultados obtidos comessa consideracao.

Para se observar o efeito da variacao de ambos os parametros nas amplitudesdos picos das tres harmonicas do espectro, fez-se uma analise buscando a distanciaao cırculo unitario conjunta dos tres polos com maiores modulos, considerando valoresmedios das simulacoes (512 casos para cada combinacao dos parametros). Criou-seum grafico tridimensional que apresenta a distancia acumulada desses polos ao cırculounitario, os quais representam as tres harmonicas procuradas, em funcao das variacoesde jitter e Ad utilizadas, Figura 5.30.

95


(a) (b)

Figura 5.22: Percentual de acertos medios do MSS para a terceira harmonica em 512simulacoes variando-se o jitter e Ad. (a) Margem de 5%. (b) Margem de 25%.

Figura 5.23: Valores obtidos para os parametros jitter e Ad utilizando a regra apre-sentada.

96


Figura 5.24: Erro medio do MSS para a primeira harmonica em 512 simulacoes com1024 combinacoes de valores de jitter e Ad.

Figura 5.25: Erro medio do MSS para a segunda harmonica em 512 simulacoes com1024 combinacoes de valores de jitter e Ad.

97


Figura 5.26: Erro medio do MSS para a terceira harmonica em 512 simulacoes com1024 combinacoes de valores de jitter e Ad.

Figura 5.27: Erro medio do MSS para a primeira harmonica com 1024 combinacoesdos parametros considerando-se somente polos com frequencias menores que 2 MHz.

98


Figura 5.28: Erro medio do MSS para a segunda harmonica com 1024 combinacoes dosparametros considerando-se somente polos com frequencias menores que 2 MHz.

Figura 5.29: Erro medio do MSS para a terceira harmonica com 1024 combinacoes dosparametros considerando-se somente polos com frequencias menores que 2 MHz.

99


Figura 5.30: Distancia acumulada dos modulos dos tres maiores polos das simulacoesao cırculo unitario em funcao das variacoes de jitter e Ad.

Analise da Ordem da Decomposicao AR

Para uma analise de qual a melhor ordem para ser utilizada na decomposicao ARdos ecos, fizeram-se algumas simulacoes utilizando valores fixos de Ad, jitter e MSS,variando-se a ordem em torno do valor de amostras usadas para caracterizar o MSS,dada a frequencia de amostragem, fs = 25 MHz, para se testar as hipoteses levantadasem [44].

Os valores usados foram: amplitude relativa da parte difusa Ad = 10%, jitter V =10%, tamanho da janela N = 512 e MSS de 1, 25 mm; 1, 5 mm; 2, 5 mm e 3, 125 mm,sendo as respectivas frequencias de MSS iguais a 0, 616 MHz; 0, 513 MHz; 0, 308 MHze 0, 246 MHz, e que correspondem ao numero de amostras usadas iguais a 40, 49, 80e 102.

Para o primeiro MSS, 1, 25 mm, variou-se o valor da ordem p de 5 a 100, emintervalos de 5 unidades, sendo que, para cada um deles, encontrava-se o erro mediodo MSS estimado para 128 simulacoes. O resultado e apresentado na Figura 5.31.

100


Figura 5.31: Erro medio do MSS de 1,25 mm, equivalente a 40 amostras, para 128simulacoes, variando-se a ordem da decomposicao AR de 5 a 100, de 5 em 5 unidades.


101


Como se pode ver, a ordem que resultou no menor erro foi 40, exatamente onumero de amostras usadas para caracterizar o MSS. Mas, em uma outra simulacao,Figura 5.32, variando-se p entre 26 e 64 de duas em duas unidades, o resultado foi umpouco diferente. Percebe-se que a ordem 40 ainda tem um erro pequeno, mas 50 e 54possuem erros um pouco mais baixos. Esses valores satisfazem a hipotese de [44] quediz que a ordem deve ser 25% maior do que o numero de amostras do MSS.

Para MSS igual a 1,5 mm, realizou-se o mesmo tipo de simulacao, mas variando-se a ordem entre 32 e 73, de uma em uma unidade, Figura 5.33. Neste caso, percebe-seque o valor esperado de 49 nao teve um erro muito baixo, apesar de 52 ter o segundomais baixo. A ordem de 67 foi a que resultou em menor erro, apesar de ser 36,7 %maior que o MSS.

Para MSS igual a 2,5 mm, aconteceu o contrario, o menor erro foi gerado poruma ordem menor que o MSS, como se pode ver na Figura 5.34, para uma variacaoda ordem de 39 a 92, de 3 em 3 unidades. Tem-se que o menor valor ocorreu comordem 60, 25 % menor que o MSS. O valor esperado, 80, obteve um erro muito alto.Realizou-se outra simulacao, variando a ordem de 25 a 225, de 25 em 25 unidades,passando alem do valor maximo de 102, N/5, apresentado em [44]. Os resultados saoapresentados na Figura 5.35, onde se pode ver que os menores erros sao para ordensproximas do MSS.

Figura 5.33: Erro medio do MSS de 1,5 mm, equivalente a 49 amostras, para 128simulacoes, variando-se a ordem da decomposicao AR de 32 a 73, de 1 em 1 unidade.

102



Figura 5.35: Erro medio do MSS de 2,5 mm, equivalente a 80 amostras, para 128 si-mulacoes, variando-se a ordem da decomposicao AR de 25 a 225, de 25 em 25 unidades.

103


Para MSS igual a 3,125 mm, fez-se uma simulacao variando a ordem entre 5 e100, de 5 em 5 unidades. O resultado e apresentado na Figura 5.36. Pode-se ver que aordem com menor erro e 80, que corresponde a um valor 21,6 % menor que o MSS, quee 102 amostras. Assim, fez-se outra simulacao, variando-se a ordem entre 70 e 127, de3 em 3 unidades, mais uma vez passando o limite de N/5, Figura 5.37, apresentandoresultados similares.

Assim, conclui-se que a ordem da regressao nao possui relacao fixa com o numerode amostras usadas para caracterizar o MSS, mas, que os valores otimos estao em tornodesse numero, variando de caso para caso.


5.3 Resultados de Sinais Coletados de Phantoms

Para a analise de sinais reais, foram utilizados phantoms cilındricos feitos de gelatina,conforme apresentado na Secao 4.1. Usaram-se 10 phantoms diferentes para se mediro MSS, descritos a seguir:

1. Gelatina de 2,9 cm de altura, com grafite, com tres fios de cobre intercaladoshorizontalmente, de forma nao coplanar, o primeiro a 4,5 mm da superfıcie, osegundo a 4,0 mm do anterior e o terceiro a 2,5 mm desse ultimo.

2. Gelatina de 3,8 cm de altura, sem grafite, com quatro fios de cobre intercaladoshorizontalmente, de forma nao coplanar, o primeiro a 6,0 mm da superfıcie, o

104



segundo a 5,0 mm do anterior, o terceiro a 4,0 mm do segundo e o quarto a 2,5mm do terceiro.

3. Gelatina de 3,7 cm de altura, sem grafite, com seis laminas de acetato intercaladashorizontalmente, a primeira a 3,5 mm da superfıcie, a segunda a 6,5 mm daanterior, a terceira a 4,5 mm, a quarta a 6,0 mm, a quinta a 7,0 mm e a sexta a6,0 mm, considerando as distancias em relacao as laminas anteriores.

4. Gelatina de 3,8 cm de altura, sem grafite, com quatro camadas finas de grafiteintercaladas horizontalmente, a primeira a 7,0 mm da superfıcie, a segunda a 6,0mm da anterior, a terceira a 6,0 mm da segunda e a quarta a 5,0 mm dessaultima.

5. Gelatina de 2,8 cm de altura, sem grafite, com tres laminas de acetato intercaladashorizontalmente, a primeira a 7,5 mm da superfıcie, a segunda a 8,0 mm daanterior e a terceira a 6,0 mm dessa ultima.

6. Gelatina de 4,4 cm de altura, sem grafite, com seis fios de cobre intercaladoshorizontalmente, de forma nao coplanar, o primeiro a 5,0 mm da superfıcie, osegundo a 6,0 mm do anterior, o terceiro a 6,0 mm, o quarto a 6,0 mm, o quintoa 6,5 mm e o sexto a 3,5 mm, considerando as distancias em relacao aos fiosanteriores.

7. Gelatina de 2,5 cm de altura, com grafite, com tres fios de cobre intercaladoshorizontalmente, de forma nao coplanar, o primeiro a 6,0 mm da superfıcie, o

105


segundo a 4,0 mm do anterior e o terceiro a 4,0 mm desse ultimo.

8. Gelatina de 2,2 cm de altura, sem grafite, com tres camadas finas de grafiteintercaladas horizontalmente, a primeira a 7,0 mm da superfıcie, a segunda a 5,5mm da anterior e a terceira a 5,0 mm dessa ultima.

9. Gelatina de 1,7 cm de altura, sem grafite, com quatro fios de cobre intercaladoshorizontalmente, de forma nao coplanar, o primeiro a 4,0 mm da superfıcie, osegundo a 2,5 mm do anterior, o terceiro a 2,0 mm do segundo e o quarto a 1,5mm do terceiro.

10. Gelatina de 1,3 cm de altura, sem grafite, com tres fios de cobre intercaladoshorizontalmente, de forma nao coplanar, o primeiro a 6,0 mm da superfıcie, osegundo a 1,5 mm do anterior e o terceiro a 0,5 mm desse ultimo.

Phantom no. 1

Este phantom possui tres fios de cobre com espacamentos de 4,5; 4 e 2,5 mm, oque corresponde a um MSS de 3,67 mm. A Figura 5.38 apresenta os sinais coletadoscom transdutor de 3,5 MHz, o sinal de referencia (em cima) e o da amostra (embaixo).Assim, calcula-se a velocidade do ultra-som na amostra, conforme [39], v = 1528, 0 m/s.Pode-se ver um eco a mais, o primeiro, devido a superfıcie da amostra, sendo geralmenteutilizado nas analises, a nao ser que nao apresente amplitude suficiente.

Como a frequencia de amostragem e muito alta, fs = 50 MHz, a analise ficoucomprometida, pois o valor de fMSS seria muito baixo, necessitando de uma ordem dadecomposicao AR bem alta. Entao, como a largura de banda do transdutor (a 3 dB)nao ultrapassa 2 MHz, a frequencia mais alta do espectro seria menor que 6,5 MHz,sendo possıvel utilizar uma taxa de amostragem de 12,5 MHz, respeitando a taxa deNyquist.

Dessa forma, para facilitar as analises, realizou-se uma subamostragem (down-sampling) para reduzir a frequencia de amostragem para fs = 12, 5 MHz [30]. Isso efeito, criando-se um novo sinal, xd[n], a partir do sinal coletado, x[n], com a seguinterelacao: xd[n] = x[4n]. A Figura 5.39 apresenta esse sinal subamostrado, ja janelado aN = 512 pontos ao redor dos ecos de interesse (em cima), e a densidade espectral depotencia da envoltoria do sinal (embaixo).

Como o estudo do efeito da ordem da decomposicao AR da secao anterior naofoi conclusivo, alguns testes serao feitos para buscar a melhor ordem para analise dosinal deste phantom. Os tempos dos picos de cada eco apresentado na Figura 5.39 (emcima) sao: t1 = 130, 32 µs, t2 = 135, 92 µs, t3 = 142, 56 µs e t4 = 145, 52 µs, assim, osintervalos entre esses ecos sao: τ1 = 5, 60 µs, τ2 = 6, 64 µs e τ3 = 2, 96 µs. Entao, comoMSS = v τ/2 e NMSS = τ fs, sendo τ o valor medio de τ e NMSS o numero de amostrasusadas para descrever o MSS, tem-se: τ = 5, 07 µs, MSS = 3, 90 mm e NMSS ≈ 64.

O valor do MSS foi bem proximo do calculado pelas distancias medidas (3,67mm). Fez-se alguns testes, variando-se a ordem da decomposicao de 32 (50% menorque NMSS) a 96 (50% maior), de 1 em 1 unidade, encontrando-se o menor erro parap = 66, proximo do esperado.

106


Figura 5.38: Sinais coletados do phantom no. 1 com transdutor de 3,5 MHz, mostram-se o sinal de referencia (em cima) e o da amostra (embaixo).

Figura 5.39: Sinal do phantom no.1, N = 512, taxa de amostragem de 12,5 MHz (emcima). Densidade espectral de potencia da envoltoria do sinal (embaixo).

107


Utilizando-se esse valor de ordem, buscou-se a frequencia dos polos com os tresmaiores modulos, conforme apresentado no Capıtulo 4. A Figura 5.40 mostra os poloscom frequencias menores que 1 MHz. Os erros encontrados para as tres primeirasharmonicas foram: 22,04; 14,72 e 19,53%, que sao valores altos, mas o jitter nesse casoe aproximadamente 27,72%, portanto os erros estao proximos dos esperados segundoas simulacoes.

Figura 5.40: Polos da decomposicao AR para o sinal do phantom no.1, apresentam-seos polos de fase positiva com frequencias menores do que 1 MHz.

Como se pode ver na Figura 5.39, o quarto eco e o que apresenta maior desviocom relacao ao espacamento medio, e tambem possui amplitude muito baixa, assim,caso ele seja desconsiderado, obtem-se erros iguais a 5,84; 3,01 e 44,38%, os quais saobem menores que os anteriores, exceto para a terceira harmonica.

Phantom no. 2

Este phantom possui quatro fios de cobre com espacamentos de 6; 5; 4 e 2,5 mm,o que corresponde a um MSS de 4,37 mm. O sinal foi coletado com transdutor de 3,5MHz e a velocidade do ultra-som na amostra e v = 1524, 6 m/s. A Figura 5.41 (emcima) apresenta o sinal ja subamostrado e janelado, N = 512, e a densidade espectralde potencia da envoltoria do sinal (embaixo).

Os tempos dos picos de cada eco apresentado na figura anterior sao: t1 =119, 92 µs, t2 = 125, 36 µs, t3 = 122, 88 µs, t4 = 138, 08 µs e t5 = 142, 08 µs, assim, osintervalos entre esses ecos sao: τ1 = 5, 44 µs, τ2 = 7, 52 µs, τ3 = 5, 20 µs e τ4 = 4, 00 µs.Entao, tem-se: τ = 5, 54 µs, MSS = 4, 20 mm, NMSS ≈ 70 e o jitter estimado e 17,87%.O valor do MSS foi bem proximo do calculado pelas distancias medidas (4,37 mm).

Fez-se alguns testes variando-se a ordem da decomposicao de 35 (50% menor queNMSS) a 105 (50% maior), de 1 em 1 unidade, encontrou-se o menor erro para p = 57,18,57% menor que o numero de amostras do MSS. A Figura 5.42 (a) mostra os polos

108



da decomposicao AR com frequencias menores que 1 MHz. O erro encontrado paraas tres primeiras harmonicas foi: 0,06; 45,20 e 40,82%.

Figura 5.42: Polos da decomposicao AR do sinal do phantom no.2, polos comfrequencias menores do que 1 MHz. (a) Ordem igual a 57. (b) Ordem igual a 68.

Pode-se ver que somente o erro da primeira harmonica e baixo, entao, buscou-seuma ordem que apresenta o menor erro do conjunto, que foi p = 68; 2,86% menor queo MSS, com erros para as tres harmonicas iguais a 11,77; 13,12 e 7,35%, os quais saobem menores que os anteriores, exceto para a primeira harmonica. A Figura 5.42 (b)mostra os polos da decomposicao AR com frequencias menores que 1 MHz para essaordem.

109


Phantom no. 3

Este phantom possui seis laminas de acetato com espacamentos de 3,5; 6,5; 4,5;6,0; 7,0 e 6,0 mm, o que corresponde a um MSS de 5,58 mm. O sinal foi coletado comtransdutor de 3,5 MHz e a velocidade do ultra-som na amostra e v = 1534, 2 m/s. Oterceiro eco apresentava amplitude muito maior que os outros, entao, para compensaressa diferenca, aplicou-se um atenuacao de 9 dB nesse terceiro eco, no trecho entre 133e 137 µs aproximadamente. A Figura 5.43 (em cima) apresenta o sinal ja subamostradoe janelado, N = 640, com o terceiro eco atenuado, e a densidade espectral de potenciada envoltoria do sinal (embaixo).

Os tempos dos picos de cada eco apresentado na figura anterior sao: t1 =118, 86 µs, t2 = 123, 94 µs, t3 = 133, 10 µs, t4 = 140, 72 µs, t5 = 149, 12 µs, t6 =156, 90 µs e t7 = 164, 48 µs, assim, os intervalos entre esses ecos sao: τ1 = 5, 08 µs,τ2 = 9, 16 µs, τ3 = 7, 62 µs, τ4 = 8, 40 µs, τ5 = 7, 78 µs e τ6 = 7, 58 µs. Entao, tem-se:τ = 7, 60 µs, MSS = 5, 80 mm, NMSS ≈ 95 e o jitter estimado e 11,16%. O valor doMSS foi proximo do calculado pelas distancias medidas (5,58 mm).

Fez-se alguns testes, variando-se a ordem da decomposicao de 48 (49,47% menorque NMSS) a 142 (49,47% maior), de 1 em 1 unidade, encontrou-se o menor erro parap = 63, 33,68% menor que o numero de amostras do MSS. A Figura 5.44 (a) mostra ospolos da decomposicao AR com frequencias menores que 1 MHz. O erro encontradopara as tres primeiras harmonicas foi: 0,21; 18,88 e 66,43%.


110


Figura 5.44: Polos da decomposicao AR para o sinal do phantom no.3, polos comfrequencias menores do que 1 MHz. (a) Ordem igual a 63. (b) Ordem igual a 96.

A ordem que apresenta o menor erro do conjunto foi p = 96, 1,05% maior queo MSS. Os erros para as tres harmonicas sao: 4,82; 0,99 e 0,04%, os quais sao bemmenores que os anteriores, exceto para a primeira harmonica, mas que tambem e baixo.A Figura 5.44 (b) mostra os polos da decomposicao AR com frequencias menores que1 MHz para essa ordem.

Phantom no. 4

Este phantom possui quatro camadas de grafite com espacamentos de 7,0; 6,0 e6,0 mm, o que corresponde a um MSS de 6,33 mm. O sinal foi coletado com transdutorde 3,5 MHz e a velocidade do ultra-som na amostra e v = 1560, 3 m/s. O primeiroeco apresentava amplitude muito maior que os outros, entao, para compensar essadiferenca, aplicou-se um atenuacao de 9 dB neste eco, no trecho entre 117 e 121 µs.O eco da quarta camada do phantom nao teve amplitude suficiente para se destacardo ruıdo, sendo desprezado, e o da terceira camada apresentou amplitude reduzida,entao, aplicou-se um ganho de 9 dB no trecho entre 142 e 146 µs. A Figura 5.45 (emcima) apresenta o sinal ja subamostrado e janelado, N = 512, com o primeiro ecoatenuado e o quarto amplificado, e a densidade espectral de potencia da envoltoria dosinal (embaixo).

Os tempos dos picos de cada eco apresentado na figura anterior sao: t1 =118, 50 µs, t2 = 126, 54 µs, t3 = 136, 30 µs e t4 = 143, 18 µs, assim, os intervalos entreesses ecos sao: τ1 = 8, 04 µs, τ2 = 9, 76 µs e τ3 = 6, 88 µs. Entao, tem-se: τ = 8, 23 µs,MSS = 6, 40 mm, NMSS ≈ 103 e o jitter estimado e 12,43%. O valor do MSS foiproximo do calculado pelas distancias medidas (6,33 mm).

Fez-se alguns testes variando-se a ordem da decomposicao de 52 (49,51% menorque NMSS) a 154 (49,51% maior), de 1 em 1 unidade, encontrou-se o menor erro parap = 64, 37,86% menor que o numero de amostras do MSS. A Figura 5.46 (a) mostra ospolos da decomposicao AR com frequencias menores que 1 MHz. O erro encontradopara as tres primeiras harmonicas foi: 0,16; 18,83 e 26,92%.

111



Figura 5.46: Polos da decomposicao AR para o sinal do phantom no.4, apresentam-seos polos de fase positiva com frequencias menores do que 1 MHz. (a) Ordem igual a64. (b) Ordem igual a 78.

112


A ordem que apresenta o menor erro do conjunto foi p = 78, 24,27% menor queo MSS. Os erros para as tres harmonicas sao: 9,47; 3,62 e 8,01%, os quais sao bemmenores que os anteriores, exceto para a primeira harmonica, mas que tambem e baixo.A Figura 5.46 (b) mostra os polos da decomposicao AR com frequencias menores que1 MHz para essa ordem.

Phantom no. 5

Este phantom possui tres folhas de acetato com espacamentos de 7,5; 8,0 e 6,0mm, o que corresponde a um MSS de 7,17 mm. O sinal foi coletado com transdutorde 3,5 MHz e a velocidade do ultra-som na amostra e v = 1535, 8 m/s. O primeiro esegundo ecos apresentavam amplitudes muito maiores que os outros dois, entao, paracompensar essas diferencas, aplicou-se uma atenuacao de 9 dB em cada um, nos trechosentre 130 e 133 µs e 139 e 143 µs respectivamente. A Figura 5.47 (em cima) apresentao sinal ja subamostrado e janelado, N = 512, com o primeiro e segundo ecos atenuados,e a densidade espectral de potencia da envoltoria do sinal (embaixo).


Os tempos dos picos de cada eco apresentado na figura anterior sao: t1 =131, 44 µs, t2 = 141, 20 µs, t3 = 151, 92 µs e t4 = 159, 44 µs, assim, os intervalos entreesses ecos sao: τ1 = 9, 76 µs, τ2 = 10, 72 µs e τ3 = 7, 52 µs. Entao, tem-se: τ = 9, 33 µs,MSS = 7, 20 mm, NMSS ≈ 117 e o jitter estimado e 12,95%. O valor do MSS foiproximo do calculado pelas distancias medidas (7,17 mm).

Fez-se alguns testes variando-se a ordem da decomposicao de 58 (50,43% menorque NMSS) a 176 (50,43% maior), de 1 em 1 unidade, encontrou-se o menor erro para

113


p = 66; 43,59% menor que o numero de amostras do MSS. A Figura 5.48 (a) mostra ospolos da decomposicao AR com frequencias menores que 1 MHz. O erro encontradopara as tres primeiras harmonicas foi: 0,09; 30,25 e 42,79%.


A ordem que apresenta o menor erro do conjunto foi p = 86, 26,50% menor queo MSS. Os erros para as tres harmonicas sao: 12,63; 6,93 e 12,06%, os quais sao bemmenores que os anteriores, exceto para a primeira harmonica, mas que tambem e baixo.A Figura 5.48 (b) mostra os polos da decomposicao AR com frequencias menores que1 MHz para essa ordem.

Phantom no. 6

Este phantom possui seis fios de cobre com espacamentos de 5,0; 6,0; 6,0; 6,0;6,5 e 3,5 mm, o que corresponde a um MSS de 5,50 mm. O sinal foi coletado comtransdutor de 3,5 MHz e a velocidade do ultra-som na amostra e v = 1567, 4 m/s.O primeiro, segundo e terceiro ecos apresentavam amplitudes muito maiores que osoutros, entao, para compensar essas diferencas, aplicou-se uma atenuacao de 9 dB emcada um, no trecho entre 110 e 127 µs. O setimo eco apresentava amplitude muitobaixa, entao, aplicou-se um ganho de 3 dB no trecho entre 155 e 158 µs. A Figura5.49 (em cima) apresenta o sinal ja subamostrado e janelado, N = 640, com os ecosatenuados e amplificados, e a densidade espectral de potencia da envoltoria do sinal(embaixo).

Os tempos dos picos de cada eco apresentado na figura anterior sao: t1 =110, 94 µs, t2 = 117, 00 µs, t3 = 125, 16 µs, t4 = 132, 68 µs, t5 = 142, 12 µs, t6 =150, 76 µs e t7 = 156, 36 µs, assim, os intervalos entre esses ecos sao: τ1 = 6, 06 µs,τ2 = 8, 16 µs, τ3 = 7, 52 µs, τ4 = 9, 44 µs, τ5 = 8, 64 µs e τ6 = 5, 60 µs. Entao, tem-se:τ = 7, 57 µs, MSS = 5, 90 mm, NMSS ≈ 95 e o jitter estimado e 15,54%. O valor doMSS foi proximo do calculado pelas distancias medidas (5,50 mm).

114



Fez-se alguns testes variando-se a ordem da decomposicao de 48 (49,47% menorque NMSS) a 142 (49,47% maior), de 1 em 1 unidade, encontrou-se o menor erro parap = 73, 23,16% menor que o numero de amostras do MSS. A Figura 5.50 (a) mostra ospolos da decomposicao AR com frequencias menores que 1 MHz. O erro encontradopara as tres primeiras harmonicas foi: 0,30; 168,63 e 126,83%.

A ordem que apresenta o menor erro do conjunto foi p = 82; 13,68% menor queo MSS. Os erros para as tres harmonicas sao: 9,70; 13,64 e 8,04%, os quais sao bemmenores que os anteriores, exceto para a primeira harmonica, mas que tambem e baixo.A Figura 5.50 (b) mostra os polos da decomposicao AR com frequencias menores que1 MHz para essa ordem.

Phantom no. 7

Este phantom possui tres fios de cobre com espacamentos de 6,0; 4,0 e 4,0 mm,o que corresponde a um MSS de 4,67 mm. O sinal foi coletado com transdutor de 3,5MHz e a velocidade do ultra-som na amostra e v = 1564, 2 m/s. O segundo e terceiroecos apresentavam amplitudes muito maiores que os outros dois, entao, para compensaressas diferencas, aplicou-se uma atenuacao de 9 dB em cada um, nos trechos entre 42 e45 µs e 48 e 51 µs respectivamente. O quarto eco apresentava amplitude muito baixa,entao, aplicou-se um ganho de 3 dB no trecho entre 53 e 56 µs. A Figura 5.51 (emcima) apresenta o sinal ja subamostrado e janelado, N = 512, com os ecos atenuadose amplificados, e a densidade espectral de potencia da envoltoria do sinal (embaixo).

115




116


Os tempos dos picos de cada eco apresentado na figura anterior sao: t1 =136, 34 µs, t2 = 142, 84 µs, t3 = 148, 76 µs e t4 = 154, 92 µs, assim, os intervalos entreesses ecos sao: τ1 = 6, 50 µs, τ2 = 5, 92 µs e τ3 = 6, 16 µs. Entao, tem-se: τ = 6, 19 µs,MSS = 4, 80 mm, NMSS ≈ 78 e o jitter estimado e 3,3%. O valor do MSS foi proximodo calculado pelas distancias medidas (4,67 mm).

Fez-se alguns testes variando-se a ordem da decomposicao de 39 (50% menorque NMSS) a 117 (50% maior), de 1 em 1 unidade, encontrou-se o menor erro parap = 88; 12,82% maior que o numero de amostras do MSS. A Figura 5.52 mostra ospolos da decomposicao AR com frequencias menores que 1 MHz. O erro encontradopara as tres primeiras harmonicas foi: 0,14; 2,50 e 2,64%. Essa tambem e a ordem queapresenta menor erro do conjunto, apesar de p = 78 ter erro de 0,28% para a segundaharmonica.

Figura 5.52: Polos da decomposicao AR para o sinal do phantom no.7, apresentam-seos polos de fase positiva com frequencias menores do que 1 MHz, ordem igual a 88.

Phantom no. 8

Este phantom possui tres camadas de grafite com espacamentos de 7,0; 5,5 e 5,0mm, o que corresponde a um MSS de 5,83 mm. O sinal foi coletado com transdutorde 3,5 MHz e a velocidade do ultra-som na amostra e v = 1525, 4 m/s. O quartoeco apresentava amplitude muito maior que os outros, entao, para compensar essasdiferencas, aplicou-se uma atenuacao de 3 dB no trecho entre 153 e 156 µs. O terceiroeco apresentava amplitude muito baixa, entao, aplicou-se um ganho de 3 dB no trechoentre 145 e 148 µs. A Figura 5.53 (em cima) apresenta o sinal ja subamostrado ejanelado, N = 512, com os ecos atenuados e amplificados, e a densidade espectral depotencia da envoltoria do sinal (embaixo).

Os tempos dos picos de cada eco apresentado na figura anterior sao: t1 =129, 48 µs, t2 = 137, 56 µs, t3 = 145, 80 µs e t4 = 154, 28 µs, assim, os intervalos entre

117



esses ecos sao: τ1 = 8, 08 µs, τ2 = 8, 24 µs e τ3 = 8, 48 µs. Entao, tem-se: τ = 8, 27 µs,MSS = 6, 30 mm, NMSS ≈ 104 e o jitter estimado e 1,72%.

Fez-se alguns testes variando-se a ordem da decomposicao de 52 (50% menor queNMSS) a 156 (50% maior), de 1 em 1 unidade, encontrou-se o menor erro para p = 107;2,88% maior que o numero de amostras do MSS. A Figura 5.54 (a) mostra os polos dadecomposicao AR com frequencias menores que 1 MHz. O erro encontrado para astres primeiras harmonicas foi: 0,0015; 1,92 e 69,52%.


Phantom no. 9

Este phantom possui quatro fios de cobre com espacamentos de 4,0; 2,5; 2,0 e 1,5mm, o que corresponde a um MSS de 2,50 mm. O sinal foi coletado com transdutorde 3,5 MHz e a velocidade do ultra-som na amostra e v = 1548, 5 m/s. O primeiro equinto ecos apresentavam amplitudes muito menores que os outros, entao, para com-pensar essas diferencas, aplicou-se um ganho de 9 dB nos trechos entre 148 e 150 µse 164 e 167 µs respectivamente. A Figura 5.55 (em cima) apresenta o sinal ja suba-mostrado e janelado, N = 512, com os ecos amplificados, e a densidade espectral depotencia da envoltoria do sinal (embaixo).

118




119


Os tempos dos picos de cada eco apresentado na figura anterior sao: t1 =148, 56 µs, t2 = 154, 48 µs, t3 = 157, 92 µs, t4 = 161, 36 µs e t5 = 165, 28 µs, assim, osintervalos entre esses ecos sao: τ1 = 5, 92 µs, τ2 = 3, 44 µs, τ3 = 3, 44 µs e τ4 = 3, 92 µs.Entao, tem-se: τ = 4, 18 µs, MSS = 3, 20 mm, NMSS ≈ 52 e o jitter estimado e 20,81%.

Fez-se alguns testes variando-se a ordem da decomposicao de 26 (50% menor queNMSS) a 78 (50% maior), de 1 em 1 unidade, encontrou-se o menor erro para p = 54;3,85% maior que o numero de amostras do MSS. A Figura 5.56 (a) mostra os polos dadecomposicao AR com frequencias menores que 1 MHz. O erro encontrado para astres primeiras harmonicas foi: 11,69; 18,15 e 16,02%.



Phantom no. 10

Este phantom possui tres fios de cobre com espacamentos de 6,0; 1,5 e 0,5 mm,o que corresponde a um MSS de 2,67 mm. O sinal foi coletado com transdutor de 3,5MHz e a velocidade do ultra-som na amostra e v = 1566, 1 m/s. O primeiro e terceiroecos apresentavam amplitudes muito menores que os outros, entao, para compensaressas diferencas, aplicou-se um ganho de 9 e 6 dB, nos trechos 152 e 155 µs e 161 e163 µs respectivamente. A Figura 5.57 (em cima) apresenta o sinal ja subamostrado ejanelado, N = 512, com os ecos amplificados, e a densidade espectral de potencia daenvoltoria do sinal (embaixo).

Os tempos dos picos de cada eco apresentado na figura anterior sao: t1 =153, 06 µs, t2 = 159, 46 µs, t3 = 162, 18 µs e t4 = 163, 46 µs, assim, os intervalos entre

120


Figura 5.57: Sinal do phantom no.10, N = 512, taxa de amostragem de 12,5 MHz(em cima). Densidade espectral de potencia da envoltoria do sinal (embaixo).

esses ecos sao: τ1 = 6, 40 µs, τ2 = 2, 72 µs e τ3 = 1, 28 µs. Entao, tem-se: τ = 3, 47 µs,MSS = 2, 70 mm, NMSS ≈ 44 e o jitter estimado e 56,41%.

Fez-se alguns testes variando-se a ordem da decomposicao de 22 (50% menor queNMSS) a 66 (50% maior), de 1 em 1 unidade, encontrou-se o menor erro para p = 38;13,64% menor que o numero de amostras do MSS. A Figura 5.58 (a) mostra os polosda decomposicao AR com frequencias menores que 1 MHz. O erro encontrado paraas tres primeiras harmonicas foi: 3,49; 33,94 e 9,70%.

A ordem que apresenta o menor erro do conjunto foi p = 45; 2,27% maior queo MSS. Os erros para as tres harmonicas sao: 6,56; 18,80 e 9,36%, os quais sao bemmenores que os anteriores, exceto para a primeira harmonica, mas que tambem e baixo.A Figura 5.58 (b) mostra os polos da decomposicao AR com frequencias menores que1 MHz para essa ordem.

121



122

Capıtulo 6

Discussoes e Conclusao

Neste capıtulo final, na primeira secao, serao apresentadas as discussoes dos resultadosdo Capıtulo 5 mostrando as vantagens e desvantagens da aplicacao do metodo dadecomposicao autorregressiva e rastreamento de polos em sinais simulados de ultra-som e sinais reais obtidos de phantoms. Na Secao 6.2, serao apresentadas as conclusoesobtidas com este trabalho, e na Secao 6.3, os proximos passos a serem seguidos nestalinha de pesquisa.

6.1 Discussao

Esta secao e dividida em duas partes: discussao dos resultados de sinais simulados ede sinais reais obtidos de phantoms. Os sinais simulados foram gerados utilizando-se o modelo apresentado no Capıtulo 2, conforme [21, 22, 35]. Os resultados foramencontrados por simulacoes de Monte Carlo e outros testes. Os phantoms foram feitosde forma a se assemelharem a tecidos moles reais, como o fıgado ou baco, conformeapresentado na Secao 4.1.

Sinais Simulados

No Capıtulo 5, apresentaram-se numerosos resultados obtidos com sinais simu-lados. Primeiramente, utilizaram-se quatro casos de teste para ilustrar o comporta-mento do metodo de analise nos sinais gerados com o modelo do meio. Para estes,alguns parametros eram mantidos fixos, Tabela 5.1, como a ordem da decomposicaoAR p = 64 tamanho da janela de dados N = 512 e MSS (espacamento medio dosespalhadores) igual a 2,5 mm.

Variavam-se dois parametros: amplitude percentual da parte difusa do sinal emrelacao a parte regular, Ad, que e uma medida da relacao sinal-ruıdo (SNR) do sinal,e a regularidade da parte regular, que e o parametro chamado jitter, uma medidapercentual da variacao das posicoes dos ecos regulares em relacao a um valor medio.Quatro conjuntos de valores foram fornecidos, Tabela 5.2, correspondendo aos casos:(1) alta regularidade e baixo ruıdo (melhor caso), (2) alta regularidade e alto ruıdo,(3) baixa regularidade e baixo ruıdo e (4) baixa regularidade e alto ruıdo (pior caso).

A Tabela 5.3 apresentou os resultados da analise de sinais simulados para es-ses quatro casos, sendo que se buscou encontrar a frequencia referente ao MSS, pela

123

Capıtulo 6 - Discussoes e Conclusao

frequencia do polo de maior modulo da decomposicao. Pode-se ver que (1) realmentee o melhor caso, pois o sinal possui pouco ruıdo e e quase periodico (baixo jitter), e(4) e o pior caso, pois e justamente o oposto, apresentando maior erro. Os casos (2) e(3) resultaram em erros intermediarios, aparentemente indistinguıveis pelo metodo.

Ainda na mesma analise, mas para uma maior compreensao dos efeitos dosparametros no metodo, buscou-se uma media dos resultados para 100 simulacoes dessesquatro casos, conforme a Tabela 5.4. Apesar dos valores medios dos erros destoaremcom os obtidos anteriormente, sendo o caso (1) o que teve maior erro medio, pode-sever que os desvios para os casos de baixa regularidade, (3) e (4), foram muito grandes,conforme a Figura 5.9. Assim, percebe-se que outras analises deveriam ser feitas paramelhores conclusoes.

Entao, partiu-se para simulacoes de Monte Carlo, as quais foram divididas emquatro tipos: primeiro, efeito da variacao da regularidade (jitter), segundo, efeito davariacao da amplitude percentual da parte difusa (Ad), SNR, terceiro, efeito de ambasas variacoes, e quarto, efeito da variacao da ordem da decomposicao autorregressiva.

No primeiro tipo, buscou-se realizar simulacoes com variacoes do jitter observandodiversos tipos de resultados, mantendo-se fixos os outros parametros como: MSS iguala 1,25 mm (fMSS = 0, 616 MHz), tamanho da janela de dados igual a 512 e ordem dadecomposicao igual a 64, alem de poucas variacoes na amplitude percentual da partedifusa.

Inicialmente, fez-se o jitter igual a 5, 15, 25, 35% e o Ad igual a 15 e 30%,obtendo-se resultados medios para a primeira harmonica considerando 512 simulacoes.Os resultados para harmonicas diferentes se referem as frequencias centrais dos tresprimeiros picos do espectro de potencia, conforme apresentado no Capıtulo 4, quecorrespondem respectivamente as frequencias dos tres polos da decomposicao AR demaior modulo.

Geraram-se histogramas com os valores de MSS encontrados. As Figuras 5.10a 5.13 apresentam os resultados para a primeira harmonica para os quatro valoresde jitter respectivamente. Pode-se ver que, na Figura 5.10 (a), a grande totalidadedos valores estava certo, e o resto com pouco desvio do esperado, em (b), para Ad =30%, o resultado foi parecido, mas apresentando alguns valores com desvios grandes.Ja na Figura 5.11 (a), ve-se menos da metade dos valores na coluna correta, varioscom desvios pequenos e um maior espalhamento. Em (b), tem-se um espalhamentolevemente maior do que (a).

Com jitter igual a 25%, os resultados ja estao muito espalhados para os doisnıveis de ruıdo, (a) 15% e (b) 30%. Na Figura 5.13, para jitter igual a 35%, osresultados espalham um pouco mais, mas ainda tem-se grande quantidade ao redor dovalor esperado, mas com erros maiores.

A Figura 5.14 apresenta um resultado consolidado, mostrando os erros de esti-mativa do MSS para as tres harmonicas para os dois valores de Ad, (a) 15% e (b)30%. Pode-se ver que os erros aumentam de uma harmonica para a outra. Para jitterbaixo, o erro para a primeira harmonica e menor do que 10%, mesmo para ruıdo alto,(b). Conforme se aumenta o jitter os erros sobem bastante, para 35%, este passa de100%. Para as outras harmonicas, um jitter de 15% ja e suficiente para um grandeerro. Quando se aumenta o ruıdo, ou seja, diminui-se a SNR (Ad = 30%), os errossofrem um leve deslocamento para cima, sendo mais acentuado quando o jitter e maior.

124

6.1 Discussao

Em seguida, analisou-se o efeito da variacao do SNR, variando-se o Ad em 10,25, 40 e 55%, com jitter de 10 e 30%. Procedendo-se da mesma forma que para o casoanterior, geraram-se os histogramas e depois um grafico consolidando os resultados. AsFiguras 5.16 a 5.19 apresentam os histogramas para a primeira harmonica, pode-se verque as ocorrencias no valor esperado sao menores que no caso analogo do jitter, apesardo valor do Ad ser menor que na Figura 5.10, o jitter e maior, e em (b) percebe-se umespalhamento bem grande.

Comparando-se os histogramas dos quatro casos de Ad, percebe-se uma dimi-nuicao das ocorrencias do valor esperado e um aumento do espalhamento numa frequenciamenor do que o ocorrido para o jitter, sendo que as partes (b) dessas figuras possuıamum espalhamento bem grande ja no primeiro caso. Esse comportamento mostra umadependencia maior dos acertos em relacao a variacoes do jitter, do que da SNR.

A Figura 5.15 apresenta os resultados consolidados para as variacoes de SNR, oserros foram maiores que para o jitter, mas, provavelmente, devido aos valores maioresdesse parametro do que dos valores de Ad. Para o caso de jitter igual a 30%, (b),percebe-se que o erro para Ad = 10% e em torno de 90%, sendo que na Figura 5.14,para Ad = 15%, para jitter igual a 5%, o erro foi menor que 10%.

Entao, para perceber qual dos dois parametros acarreta maiores efeitos no metodo,variou-se os dois. Utilizaram-se cinco valores de cada um, 5, 15, 25, 35 e 45%, to-talizando 25 combinacoes, para as quais foram feitas 512 simulacoes, buscando-se aquantidade de acertos para dois casos, acerto com margem de 5% e 25% ao redor dovalor esperado, conforme [22, 35].

Assim, foi obtida a Figura 5.20, para a primeira harmonica, que apresenta duaspartes, (a), margem de 5% e, (b), margem de 25%, sendo que cada uma, possui cincocurvas, uma para cada valor de Ad, variando-se no eixo das abscissas, o jitter, e no eixodas ordenadas, o percentual de acertos. Pode-se ver que, para margem de 5%, quantomaior o Ad (menor SNR), as curvas apresentam menores valores de acertos para todosos jitter, sendo que para Ad igual a 5 e 15%, percebe-se uma variacao maior quando sevaria o jitter, principalmente no primeiro caso.

Para margem de 25%, ocorre um deslocamento do percentual de acertos paracima, como era de se esperar, mas percebe-se uma sensibilidade maior a variacoesdo jitter. Com essa figura, percebe-se que o espalhamento dos resultados aumentaconforme se aumenta o jitter e a SNR, mas nao fica muito claro a qual o sistematem maior sensibilidade. As Figuras 5.21 e 5.22 mostram resultados analogos para asegunda e terceira harmonicas respectivamente, as quais sao muito mais sensıveis avariacoes de jitter e SNR, apresentando queda do percentual de acertos muito maisacentuada do que a primeira, quando se varia qualquer um dos parametros.

Para uma melhor visualizacao destes efeitos, fez-se outra simulacao, variando-seos parametros em 32 valores cada, resultando em 1024 combinacoes. Essas variacoesforam nao-lineares, como mostra a Equacao 5.1, variando o jitter e o Ad de 1 a 38%.Para cada uma, fez-se 512 simulacoes, apresentado-se os resultados de erros medios deestimativa nas Figuras 5.24 a 5.26 para cada uma das harmonicas. Nessas, fica claroque e ao jitter que o sistema e mais sensıvel, pois, para a primeira harmonica, a partirde 8% de jitter, tem-se um aumento grande do erro, mas e so a partir de 24% de Adque o mesmo ocorre. Resultados analogos sao apresentados para as outras harmonicas,que resultam em erros maiores para os mesmos casos citados.

125


Foi feito um pequeno ajuste no algoritmo que acarretou uma melhora significativanos resultados. Como o maior valor de fMSS procurado era 1,848 MHz para a terceiraharmonica, pois o MSS era 1,25 mm, nao se fazia necessario buscar valores numafaixa acima de 2 MHz. Dessa forma, limitando-se a area de procura, obtiveram-se osresultados apresentados nas Figuras 5.27 a 5.29, percebendo-se uma diminuicao grandedo erro que caiu para um maximo de aproximadamente 30% para a primeira harmonicaem relacao a 150% obtido anteriormente. Para as outras harmonicas, obtiveram-seresultados ainda melhores, 20% em relacao a 200%, e 15% em relacao a 250%, para asegunda e terceira harmonicas respectivamente.

Um outro resultado interessante e a variacao dos modulos dos tres maiores polosda decomposicao que correspondem as tres harmonicas, com respeito a mesma variacaodos parametros, Figura 5.30. Pode-se ver que, conforme previsto teoricamente em [22],as amplitudes dos picos do espectro diminuem conforme se aumenta o jitter e o ruıdo,pois as distancias do cırculo unitario aumentam. Assim, percebe-se que, com o aumentodo ruıdo e/ou diminuicao da regularidade (aumento do jitter), a energia espectral seespalha, diminuındo a altura dos picos, deixando-os mais distribuıdos na frequencia.

Em seguida, realizaram-se testes com variacoes da ordem da decomposicao AR.Utilizaram-se alguns valores de MSS e variou-se a ordem, mantidos constantes osparametros jitter e Ad, iguais a 10%, alem do tamanho da janela de dados, 512 pontos.Os valores de MSS usados foram 1,25; 1,5; 2,5 e 3,125 mm.

Para o primeiro valor, 1,25 mm, variou-se a ordem de 5 a 100, de 5 em 5 unidades,buscando-se qual delas apresentava um erro menor de estimativa de MSS medio em128 simulacoes. O resultado, Figura 5.31, foi 40, que equivalia ao numero de amostrasusado para caracterizar o MSS, considerando a frequencia de amostragem igual a 25MHz. Na Figura 5.32, pode-se ver que os valores de ordem iguais a 50 e 54 tambemapresentaram erros bem proximos, que sao valores condizentes com o resultado de [44],o qual diz que a ordem deve ser 25% maior que o MSS.

Para MSS = 1, 5 mm, o resultado para simulacoes com ordem entre 32 a 73, de 1em 1, Figura 5.33, foi 67; 36,7 % maior que o valor do MSS, 49 amostras. Mas, paraa ordem 52, obteve-se o segundo menor erro. Para MSS = 2, 5 mm, o resultado ja foidiferente, pois, variando-se a ordem entre 39 a 92, de 3 em 3, Figura 5.34, o menor errofoi para ordem 60, que corresponde a um numero de amostras 25% menor que o doMSS usado nas simulacoes, que seria 80, para o qual obteve-se um erro bem grande.

Para MSS = 3, 125 mm, fez-se duas simulacoes, uma com ordem entre 5 e 100, de5 e 5, e outra variando de 70 a 127, de 3 em 3, o menor valor de erro foi para ordem 80;21,6% menor que o numero de amostras usadas para caracterizar o MSS, 102. Segundo[44], nao se deve usar ordem maior que N/5, sendo N o tamanho da janela, o queacarretaria erros elevados, como foi o caso da Figura 5.35, para MSS igual a 2,5 mm,mas para valores um pouco acima de N/5, em torno de N/4, assim, deve-se usar umtamanho de janela condizente com o MSS que se deseja detectar.

Os resultados para variacao da ordem da decomposicao AR nao foram conclusivos,pois nao se encontrou uma regra clara para a determinacao desse parametro, como foiobtido em [44]. Dessa forma, para as analises dos sinais reais de phantoms, fez-se algunstestes antes de se escolher a ordem a ser usada. Foi possıvel observar que os valoresotimos de ordem possuem variacao maxima de entre 50% menor e 50% maior do queo nuumero de amostras usadas para caracterizar o MSS.

126

6.1 Discussao

Sinais Reais

Foram analisados sinais de dez phantoms feitos com gelatina, conforme apresen-tado na Secao 5.3. Nessas analises, buscaram-se os polos que possuıssem frequenciasinferiores a 1 MHz, pois os valores de MSS usados eram tais que a terceira harmonicatinha frequencia menor que esse valor. Tambem, ao se encontrar os tres polos commaior modulo, esses foram organizados conforme a frequencia, pois, devido a grandequantidade de ruıdo presente nos sinais e ao alto jitter de alguns deles, nem sempre ospolos possuıam modulos decrescentes, porque o espectro podia ter um pico maior parauma frequencia mais elevada, como se pode ver nas figuras da Secao 5.3.

Utilizou-se, tambem, uma subamostragem para diminuir a frequencia de amos-tragem de 50 MHz para 12,5, pois o transdutor usado possuıa frequencia central de 3,5MHz e largura de banda em torno de 2 MHz, assim a maior frequencia do espectrodos ecos seria menor que 6,5 MHz, portanto e preferıvel trabalhar com uma frequenciamenor. O transdutor de 8 MHz tambem foi usado, mas os sinais apresentaram menorqualidade, provavelmente, porque, nessa frequencia, os meios usados como interfacenao apresentavam ecos com amplitudes suficientes para serem analisados.

As janelas de dados usadas foram de 512 pontos, exceto para os phantoms numero3 e 6, sendo usado 640, pois estes possuıam sete ecos, precisando de uma janela maiorpara englobar todos eles. Utilizou-se, tambem, uma amplificacao ou atenuacao emalguns ecos de algumas amostras, para que todos possuıssem amplitudes semelhantes,facilitando as analises.

Faziam-se simulacoes usando varios valores de ordem da decomposicao AR, va-riando entre 50% menor e maior que o numero de amostras usado para caracterizaro MSS nos sinais, buscando-se o valor que obtivesse o menor erro de estimativa paraa primeira harmonica. Entao, encontrava-se os erros para as tres harmonicas comessa ordem. Buscou-se tambem o valor que resultava no menor erro do conjunto dasharmonicas.

O primeiro phantom possuıa grafite em sua constituicao e tres fios de cobre nomeio. O MSS e a velocidade do ultra-som estimados pela inspecao visual do sinal deeco foram 3,90 mm e 1528,0 m/s respectivamente, assim, tem-se fMSS = 0, 196 MHz.O numero de amostras usadas para caracterizar o MSS no sinal foi 64.

Neste sinal, nao foram feitas compensacoes, assim, o quarto eco apresentava umaamplitude relativamente baixa, sendo tambem bem proximo do terceiro, resultando emalto jitter do sinal. Dessa forma, fez-se duas analises, a primeira considerando esse ecoe a segunda, nao. Para a primeira, com MSS = 3, 90 mm, os erros obtidos, usando aordem igual a 66, foram 22,04; 14,72 e 19,53%, que sao valores condizentes, segundoas simulacoes com sinais artificiais, para jitter em torno de 28%.

Para a segunda analise, o MSS estimado, sem o quarto eco, foi 4,68 mm, obtendo-se erros iguais a 5,84; 3,01 e 44,38 %, claramente valores bem mais baixos do que osanteriores, exceto para a terceira harmonica, a qual e mais sujeita a variacoes no sinal,nem estando aparente na Figura 5.39, a densidade espectral de potencia do eco.

A partir do segundo phantom, comecou-se a utilizar a equalizacao dos ecos, con-forme apresentado na Secao 5.3. Nesse phantom, feito sem grafite e com quatro fios decobre, o MSS e velocidade estimados no meio pela inspecao visual foram 4,20 mm e1524,6 m/s respectivamente, assim, tem-se fMSS = 0, 181 MHz. O numero de amos-

127


tras usadas para caracterizar o MSS no sinal foi 70, sendo o jitter em torno de 18%.Os erros encontrados para a ordem de 57, que era a que produzia o menor erro daprimeira harmonica, foram 0,06; 45,20 e 40,82%.

Pode-se ver que o primeiro valor e muito baixo, bem menor que 1%, mas os errospara as outras harmonicas estao altos. Assim, utilizou-se a ordem 68, encontrando-sevalores de erro iguais a 11,77; 13,12 e 7,35%, os quais sao valores bons para o nıvelde jitter. Pode-se ver que, escolhendo-se a ordem certa, os erros obtidos sao baixos,ate melhores que os obtidos na simulacao com ecos artificiais, principalmente para asegunda e terceira harmonicas. Pode-se ver, Figura 5.42, que os polos para o primeirocaso, (a), nao se comportam como o esperado, ja em (b), tem-se o que se espera, umdecaımento quase ideal (quase, pois o segundo polo tem modulo um pouco menor, maso algoritmo ainda o considera como referente a segunda harmonica, como explicadoanteriormente).

Para o terceiro phantom, feito sem grafite e com seis laminas de acetato, o MSSe a velocidade estimados foram 5,80 mm e 1534,2 m/s respectivamente, assim, tem-sefMSS = 0, 132 MHz. O numero de amostras usadas para caracterizar o MSS no sinalfoi 95, sendo o jitter em torno de 11%. Os erros encontrados para ordem 63 foram0,21; 18,88 e 66,43%.

Mais uma vez, pode-se ver que o primeiro valor e muito baixo, menor que 1%,mas os erros para as outras harmonicas estao altos. Assim, utilizou-se a ordem 95,encontrando-se valores de erro iguais a 4,82; 0,99 e 0,04%, os quais sao valores muitobons para o nıvel de jitter. Neste caso, os erros para as outras harmonicas foram bempequenos, pode-se ver pela Figura 5.44 (b) que os polos se comportam praticamentecomo o esperado. Isso e comprovado pela Figura 5.43, que mostra o PSD com os tresprimeiros picos mais proeminentes.

Para o quarto phantom, com quatro camadas finas de grafite como interfaces, oMSS e a velocidade estimados foram 6,40 mm e 1560,3 m/s respectivamente, assim,tem-se fMSS = 0, 122 MHz. O numero de amostras usadas para caracterizar o MSSno sinal foi 103, sendo o jitter em torno de 12,5%. Os erros encontrados para ordem64 foram 0,16; 18,83 e 26,92%.

Para se obter o menor erro do conjunto de harmonicas, utilizou-se a ordem 78,encontrando-se valores de erro iguais a 9,47; 3,62 e 8,01%, os quais sao valores bonspara o nıvel de jitter, apesar de serem maiores que os do phantom no. 2. Os polosda decomposicao para as duas ordens comportaram-se conforme o esperado, como sepode ver pela Figura 5.46.

Para o quinto phantom, feito sem grafite e com tres laminas de acetato, o MSSe a velocidade estimados foram 7,20 mm e 1535,8 m/s respectivamente, assim, tem-sefMSS = 0, 107 MHz. O numero de amostras usadas para caracterizar o MSS no sinalfoi 117, sendo o jitter em torno de 13%. Os erros encontrados para ordem 66 foram0,09; 30,25 e 42,79%.

Para se obter o menor erro do conjunto de harmonicas, utilizou-se a ordem 86,encontrando-se valores de erro iguais a 12,63; 6,93 e 12,06%, os quais sao valores bonspara o nıvel de jitter, apesar de serem maiores que os do phantom no. 2. Os polos dadecomposicao para as duas ordens comportaram-se conforme o esperado, da mesmaforma que para o phantom anterior, como se pode ver pela Figura 5.48.

Para o sexto phantom, feito sem grafite e com seis fios de cobre, o MSS e a

128

6.1 Discussao

velocidade estimados foram 5,90 mm e 1567,4 m/s respectivamente, assim, tem-sefMSS = 0, 133 MHz. O numero de amostras usadas para caracterizar o MSS no sinalfoi 95, sendo o jitter em torno de 15,5%. Os erros encontrados para ordem 73 foram0,30; 168,63 e 126,83%. Neste caso, para as outras harmonicas, os erros foram muitograndes, mas que o dobro do esperado.

Para se obter o menor erro do conjunto de harmonicas, utilizou-se a ordem 82,encontrando-se valores de erro iguais a 9,70; 13,64 e 8,04%, os quais sao valores bonspara o nıvel de jitter, e bem menores para as outras harmonicas do que os anteriores,percebe-se que neste caso, a ordem influe muito nos resultados. Apesar do PSD ter umbom comportamento, Figura 5.49, percebem-se picos altos depois do terceiro, e pelospolos da decomposicao para a primeira ordem, Figura 5.50 (a), entende-se os errosaltos, devido a polos de altas frequencias com modulo grandes. Mas, para a segundaordem, o comportamento destes foi quase como esperado, Figura 5.50 (b).

Para o setimo phantom, feito com grafite e com tres fios de cobre, o MSS e avelocidade estimados foram 4,80 mm e 1564,2 m/s respectivamente, assim, tem-sefMSS = 0, 163 MHz. O numero de amostras usadas para caracterizar o MSS no sinalfoi 78, sendo o jitter em torno de 3,3%, que e um valor baixo. Os erros encontradospara ordem 88 foram 0,14; 2,50 e 2,64%, os quais sao resultados muito bons, condi-zentes com o jitter, conforme simulado para sinais artificiais, mostrando a dependenciadesse parametro como concluıdo anteriormente. Essa ordem tambem e a que resultaem menor erro de conjunto. A Figura 5.51 mostra que o PSD, que possuı os tres pri-meiros picos bem definidos e com amplitudes progressivamente menores e com largurasmaiores. A Figura 5.52 mostra os polos da decomposicao exatamente como esperado.

Para o oitavo phantom, com tres camadas finas de grafite, o MSS e a veloci-dade estimados foram 6,30 mm e 1525,4 m/s respectivamente, assim, tem-se fMSS =0, 121 MHz. O numero de amostras usadas para caracterizar o MSS no sinal foi 104,sendo o jitter em torno de 1,7%. Os erros encontrados para ordem 107 foram 0,0015;1,92 e 69,52%. Este e o menor erro de primeira harmonica obtido para os sinais de phan-tom, o jitter e bem baixo, como esperado para um erro desses. A terceira harmonicaapresentou erro alto.

Para se obter o menor erro do conjunto de harmonicas, utilizou-se a ordem 103,encontrando-se valores de erro iguais a 0,19; 0,75 e 0,32%, os quais sao valores bons,condizentes com o nıvel de jitter. A terceira harmonica possui erro bem menor que oanterior. O PSD tem um bom comportamento, Figura 5.53, os polos da decomposicaopara a primeira ordem, Figura 5.54 (a), nao se comportam muito bem, mas em (b)tem-se o esperado.

Para o nono phantom, feito sem grafite e com quatro fios de cobre, o MSS ea velocidade estimados foram 3,20 mm e 1548,5 m/s respectivamente, assim, tem-sefMSS = 0, 242 MHz. O numero de amostras usadas para caracterizar o MSS no sinalfoi 52, sendo o jitter em torno de 21%. Os erros encontrados para ordem 54 foram11,69; 18,15 e 16,02%, os quais sao baixos, quando comparados com os valores obtidospara esse nıvel de jitter nas simulacoes.

Para se obter o menor erro do conjunto de harmonicas, utilizou-se a ordem 45,encontrando-se valores de erro iguais a 11,74; 7,58 e 8,74%, os quais sao valores bons,condizentes com o nıvel de jitter, sendo metade dos obtidos anteriormente para asoutras harmonicas. Pode-se ver na Figura 5.55, que o PSD nao tem um bom com-

129


portamento. Devido ao jitter maior, a largura dos picos com maior frequencia saomaiores, resultando em um cruzamento entre estes, prejudicando a analise, conformeprevisto teoricamente em [22]. Os polos da Figura 5.56 (b) apresentam um bom com-portamento, sendo que aqueles correspondendo as outras harmonicas tem modulos bemmenores que a primeira, como era de se esperar para o nıvel de jitter.

Para o decimo phantom, feito sem grafite e com tres fios de cobre, o MSS e avelocidade estimados foram 2,70 mm e 1566,1 m/s respectivamente, assim, tem-sefMSS = 0, 290 MHz. O numero de amostras usadas para caracterizar o MSS no sinalfoi 44, sendo o jitter em torno de 56,5%, o qual e um valor muito alto. Os errosencontrados para ordem 38 foram 3,49; 33,94 e 9,70%, os quais sao valores muito bons,quando comparados com os valores obtidos para esse nıvel de jitter nas simulacoes.Percebe-se que, mesmo com o ruıdo alto, como o jitter era baixo, os resultados forambons.

Para se obter o menor erro do conjunto de harmonicas, utilizou-se a ordem 45,encontrando-se valores de erro iguais a 6,56; 18,80 e 9,36%, os quais sao, mais umavez, valores bons, quando se consideram as simulacoes para esse nıvel de jitter. Assim,percebe-se que, com um bom controle da ordem, pode-se encontrar erros de estimativamenores do que os apresentados nas simulacoes.

Pode-se ver na Figura 5.57, que o PSD, como para o phantom anterior, nao temum bom comportamento, devido ao jitter maior, sendo maiores as larguras dos picoscom maior frequencia, resultando em um cruzamento entre estes. Os polos dos doiscasos da Figura 5.57 se comportam muito mal, o que ocasiona erros mais elevados, masdevido a organizacao dos polos detectados segundo a frequencia, esses problemas saoatenuados. Tambem, neste phantom e no anterior, como o MSS era menor que nosoutros, o numero de polos com frequencias menores que 1 MHz ficou reduzido a tresou quatro, como na Figura 5.56 (a), o que mascara os resultados.

Como se pode ver nas analises de sinais de phantoms, a escolha da ordem daregressao AR foi de extrema importancia para um bom resultado. Pode-se ver que ovalor que resultava no menor erro estava proximo do numero de amostras usadas paracaracterizar o MSS, mas nao era regido por uma regra clara, portanto so podia serobtido por testes.

6.2 Conclusao

Neste trabalho, foram realizados estudos de sinais de ultra-som, tanto artificiais, gera-dos por modelos simplificados de tecidos biologicos, quanto reais, obtidos por inspecoesde ultra-som em phantoms feitos de gelatina, grafite, fios de cobre e laminas de acetato.Os modelos utilizados consideravam que os meios, tecidos moles como fıgado ou baco,podiam ser simplificados por linhas de espalhadores de pulsos de ultra-som. Essaslinhas eram compostas de dois tipos de espalhadores, um aleatorio, que modelava aparte difusa do tecido, e outra regular, com distancias definidas por uma propriedadedo meio chamada espacamento medio dos espalhadores (MSS).

Esses espalhadores nem sempre eram regularmente espacados o que era caracteri-zado pelo jitter do meio, ou seja, o quanto esse espacamento divergia da media. A partedifusa tambem podia apresentar variacoes, caracterizadas no modelo pela amplitude

130

6.2 Conclusao

da parte difusa do eco retroespalhado em relacao a parte regular.

Os sinais foram analisados usando decomposicao autorregressiva e rastreamentodos polos. O princıpo da analise era que a propriedade de MSS podia ser caracterizadapelos picos do espectro de potencia do eco ao quadrado, entao, os polos da decomposicaoestavam diretamente relacionados a esses picos. Dessa forma, era possıvel encontraro MSS sem a necessidade de estimar o PSD, nem enfrentar os problemas decorrentesdisso. O trabalho se reduzia a utilizar ferramentas conhecidas e bem estabelecidasna literatura, como LMS ou metodo da autocorrelacao, para encontrar os coeficientesda decomposicao AR. Em seguida, utilizava-se um metodo de busca de raızes de po-linomios, como a solucao pelos autovalores da matriz companheira, para encontrar ospolos.

De posse dos polos, utilizavam-se tecnicas de rastreamento para encontrar aquelesque possuem os modulos mais proximos do cırculo unitario, pois, estes sao os queproduzem os maiores picos do espectro de potencia. Encontrando-se a frequencia dessespolos, tinha-se a frequencia central do pico do espectro. Usaram-se os tres primeirospicos do PSD, correspondendo a harmonicas de fMSS, para encontrar o valor de MSS.

Para validar o metodo, varias simulacoes foram feitas utilizando sinais artificiaispelo metodo de Monte Carlo. Fez-se testes para se descobrir o efeito de variacoesdos parametros: jitter, amplitude da parte difusa (Ad) e ordem da decomposicao AR.Percebeu-se que o aumento do jitter ou do Ad provocava um aumento do erro deestimativa do MSS, mas o metodo e mais sensıvel ao primeiro parametro, o qual apartir de 8% ja resultava em uma degradacao dos resultados.

A respeito da ordem da decomposicao AR, nao foram obtidos resultados con-clusivos, descobriu-se, como se observa na literatura, [44, 22], uma grande influenciadesse parametro nos resultados, mas nao se encontrou uma regra exata para a esco-lha. Percebeu-se que valores entre 50% menores ou maiores que o numero de amostrasusadas para caracterizar o MSS nos sinais apresentavam resultados satisfatorios. En-tretanto, chegou-se a conclusao que testes deveriam ser feitos com cada sinal a seranalisado para se encontrar a ordem adequada, sendo este o procedimento adotado nasanalises de sinais reais.

Para a analise dos sinais reais, utilizaram-se dez phantoms, sendo que para cadaum, fez-se simulacoes para buscar a ordem adequada da decomposicao AR, ou seja,aquela que produzia o menor erro de estimativa do MSS. Utilizaram-se duas abor-dagens, primeira, a ordem que minimizava o erro da primeira harmonica e, segunda,a que minimizava o erro do conjunto. Fez-se necessario tambem o uso de um pre-processamento dos sinais, para equalizar as amplitudes dos ecos, que, as vezes, erammuito maiores ou menores que os outros.

Os sinais apresentavam nıveis de ruıdo diversos, mas esse nao foi quantificado, poisseria um procedimento bastante complicado. Em relacao ao jitter, que e o parametromais sensıvel do metodo, foram obtidos varios valores, os quais cobriam os casos apre-sentados nas simulacoes com sinais artificiais.

Dessa forma, os resultados foram muito bons, pois os erros chegavam a ser meno-res que os apresentados nas simulacoes, e, na maioria dos casos, eram compatıveis comos resultados obtidos anteriormente. Assim, pode-se dizer que o modelo usado para ostecidos, apesar de simples, foi fiel ao encontrado nos phantoms, da mesma forma queem [22], o qual tambem apresentou bons resultados para sinais de fıgado bovino.

131


Pode-se concluir, tambem, que o metodo da decomposicao autorregressiva foieficiente na deteccao do MSS dos sinais, comportando-se tao bem quanto os metodosutilizados em [35], que usou a magnitude e fase do espectro de potencia para estimaro MSS, em [45], que usou a tecnica da autocorrelacao espectral, e em [22], que usou aanalise de espectro singular.

A grande vantagem deste metodo e ser simples, nao precisar de estimativa doespectro de potencias do eco e encontrando resultados similares a esses metodos citados.Ele tambem traz uma contribuicao aos estudos de ecos de ultra-som retroespalhadospor meios, que e a relacao entre as amplitude dos primeiros picos do PSD, que segundodesenvolvido teoricamente em [22], devem apresentar amplitudes menores e largurasde banda maiores conforme se aumenta o jitter e/ou o a amplitude da parte difusa domeio, o que foi observado para o caso do jitter nas analises dos sinais de phantoms.

6.3 Proximos Passos

Os proximos passos deste trabalho seriam analises mais aprofundadas de sinais dephantoms utilizando ecos obtidos em varreduras desses meios, pois assim, seria possıvelfazer analises estatısticas dos resultados, da mesma forma que foi feito pelo metodo deMonte Carlo para os sinais artificiais.

Seria interessante, tambem, testar alguns metodos mais refinados de pre-processamentodos sinais reais, buscando eliminar alguns problemas encontrados, como a diferenca dosecos. Utilizar outras formas de construcao dos phantoms, usando outros materiais, etentando obter um controle maior do espacamento produzido. Os resultados seriammelhorados se fosse possıvel usar um sistema mais preciso de inspecao por ultra-som,com transdutores de alta qualidade.

Faz-se necessario, nos proximos trabalhos, analisar melhor os efeitos da ordem dadecomposicao AR nos resultados. Buscando-se uma regra mais precisa para a escolhadesse parametro, ja que se pode ver a sua importancia para uma boa estimativa doMSS.

Poderia-se aumentar os conhecimentos do comportamento do metodo, pela aplicacaodeste a sinais obtidos de tecidos moles reais, como amostras de fıgado ou baco, de formaa se comparar os resultados obtidos com phantoms, podendo validar definitivamente ometodo e o modelo usados.

Alguns outros metodos apresentados na literatura, como o uso da entropia citadopor [22] ou a filtragem de kalman aplicada a series temporais [40], poderiam ser usadosde forma a se aprofundar ainda mais o conhecimento desta area. Tambem, buscando-serelacionar os resultados de estimativa de MSS de amostras normais e patologicas.

132

Referencias Bibliograficas

[1] ACHENBACH, Jan D. Modeling for Quantitative Non-Destructive Evaluation.Ultrasonics, v.40, p.1-10, 2002.

[2] BENCHEQROUNE, A.; BENSEDDIK, M; HAJJARI, A. Tracking of Time-Varying Frequency of Sinusoidal Signals. Signal Processing, v.78, p.191-199, 1999.

[3] BRIDAL, S. L.; SAIED, A.; CHERIN, E.; LEFEBVRE, F.; LAUGIER, P.; BER-GER, G. High-Resolution Quantitative Echography with Tissue Parameters. Jour-nal d’Echographie et de Medecine par Ultrasons, v.19, n.2/3, 1998, p.204-211.

[4] COHEN, Fernand S.; GEORGIOU, Georgia; HALPERN, Ethan J. WOLD De-composition of the Backscatter Echo in Ultrasound Images of Soft Tissue Organs.IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control, v.44,n.2, p.460-471, 1997.

[5] DINIZ, Paulo R. S. Adaptive Filtering: Algorithms and Practical Implementation.The Kluwer International Series in Engineering and Computer Science. KluwerAcademic Publishers, Boston, 1999.

[6] EL-JAROUDI, Amro; MAKHOUL, John. Speech Analysis Using Discrete Spec-tral Modeling. In IEEE Proceedings of the 32nd Midwest Symposium on Cir-cuits and Systems, 1989. Disponıvel em http://intl.ieeexplore.ieee.org/

xpl/abs_free.jsp?arNumber=101800, acessado em 09/12/2004.

[7] ENDERLE, John, BLANCHARD, Susan, BRONZINO, Joseph. Introduction toBiomedical Engineering. Academic Press, New York, 2000. Cap.15: Ultrasound.

[8] ERKAMP, Ramon Q.; SKOVORODA, Andrei R.; et al. Measuring the NonlinearElastic Properties of Tissue-like Phantoms. IEEE Transactions on Ultrasonics,Ferroelectrics, and Frequency Control, v.51, n.4, p.410-419, 2004.

[9] FISH, P. Physics and Instrumentation of Diagnostical Medical Ultrasound. JohnWiley & Sons, New York, 1990.

[10] GALASSI, Mark; DAVIES, Jim; THEILER, James; GOUGH, Brian; JUNGMAN,Gerard; BOOTH, Michael; ROSSI, Fabrice. GNU Scientific Library - ReferenceManual. 1.5 ed. Free Software Foundation, 2004.

[11] GEDDES, L. A., BAKER, L.E. Principles of Applied Biomedical Instrumentation.3rd ed. John Wiley & Sons, Wiley-Interscience, New York, 1989.

133

REFERENCIAS BIBLIOGRAFICAS

[12] GOLUB, Gene H.; LOAN, Charles F. Van. Matrix Computations. 3rd ed. JohnsHopkins University Press, Johns Hopkins Series in the Mathematical Sciences,Baltimore, 1996.

[13] HAYES, Monson H. Statistical Digital Signal Processing and Modeling. JohnWiley & Sons, Inc., New York, 1996.

[14] HAYKIN, Simon. Adaptive Filter Theory. 4th ed. Prentice Hall Information andSystem Science Series. Prentice Hall, Inc., New Jersey, 2001.

[15] HUISMAN, Henkjan J.; THIJSSEN, Johan M. An In Vivo Ulrasonic Model of Li-ver Parenchyma. IEEE Transactions on Ultrasonics, Ferroelectrics, and FrequencyControl, v.45, n.3, p.739-750, 1998.

[16] ISHIMARU, A. Wave Propagation and Scattering in Random Media. AcademicPress, Londres, 1978.

[17] JENDERKA, Klaus V.; GAERTNER Tilo; et al. System Independent TissueTyping of Human Testis and Prostate. Proc. IEEE Intern. Ultrasonics Symposium,Lake Tahoe, Nevada USA, 1999.

[18] KADAH, Y. M., FARAG, A. A. et al. Classification Algorithms for QuantitativeTissue Characterization of Diffuse Liver Disease from Ultrasound Images. IEEETransactions on Medical Imaging, v. 15, n. 4, p. 466-478, 1996.

[19] KAY, Steven M. Modern Spectral Estimation: Theory & Application. PrenticeHall Signal Processing Series. Prentice Hall, Inc., New Jersey, 1988.

[20] KINSLER, Lawrence E.; FREY, Austin R.; COPPENS, Alan B.; SANDERS,James V. Fundamentals of Acoustics. 3rded. John Wiley & Sons, Inc., New York,1984.

[21] LANDINI, Luigi, VERRAZZANI, Lucio. Spectral Characterization of TissuesMicrostructure by Ultrasounds: A Stochastic Approach. IEEE Transactions onUltrasonics, Ferroelectrics and Frequency Control, v.37, n. 5, p. 448-456, 1990.

[22] MACIEL, Carlos Dias. Analise de Espectro Singular Aplicada a Sinais Ultra-Sonicos. Tese (Tese de Doutorado em Ciencias em Engenharia Biomedica). Ori-entador: Wagner Coelho de Albuquerque Pereira. Universidade Federal do Rio deJaneiro, Rio de Janeiro, 2000.

[23] MACIEL, Carlos Dias; PEREIRA, Wagner Coelho de Albuquerque. Modelagem eProcessamento e Speckle em Imagem Ultra-Sonica: Uma Revisao. RBE - Cadernode Engenharia Biomedica, v.13, n.1, p.71-89, 1997.

[24] MAINARDI, Luca T.; BIANCHI, Anna M.; BASELLI, Giuseppe; CERUTTI,Sergio. Pole-Tracking Algorithms for the Extraction of Time-Variant Heart RateVariability Spectral Parameters. IEEE Transactions on Biomedical Engineering,v.42, n.3, p.250-259, 1995.

134


[25] MAINARDI, L. T.; PETRUCCI, E.; BIANCHI, A. M.; BALIAN, V.; BERTI-NELLI, M.; MAINARDI, M.; CERUTTI, S. Time-Variant Spectral Analysis forthe Studying of Transient Ischemia during Dipyridamole Stress Test. In IEEEConf. Computers in Cardiology, 1994. Disponıvel em http://intl.ieeexplore.

ieee.org/xpl/abs_free.jsp?arNumber=470261, acessado em 09/12/2004.

[26] MARON, Melvin J. Numerical Analysis: A Practical Approach. 2nd ed. CollierMacmillan, 1987.

[27] MARPLE, S. Lawrence, Jr. Computing the Discrete-Time Analytic Signal viaFFT. IEEE Transactions on Signal Processing, v.47, n.9, p.2600-2603, 1999.

[28] NAIDU, V. P. S.; REDDY, M. R. S. Time-Variant Power Spectral Analysis ofHeart-Rate Time Series by Autoregressive Moving Average (ARMA) Method. TheIndian Academy of Sciences, Bangalore. Academy Proceedings in Engineering Sci-ences, Sadhana, v.28, part.6, p.1027-1035, 2003.

[29] NARAYANAN, V. Manoj; MOLTHEN, Robert C.; et al. Studies on UltrasonicScattering from Quasi-Periodic Structures. IEEE Transactions on Ultrasonics,Ferroelectrics and Frequency Control, v.44, n.1, p. 114-124, 1997.

[30] OPPENHEIM, Alan V.; SCHAFER, Ronald W.; BUCK, John R. Discrete-TimeSignal Processing. Prentice Hall Signal Processing Series. Prentice-Hall, Inc., NewJersey, 1999.

[31] PAPOULIS, Athanasios. Probability, Random Variables, and Stochastic Processes.3rd ed. McGraw-Hill, Inc., Singapore, 1991.

[32] PATOMAKI, L.; KAIPIO, J. P.; KARJALAINEN, P. A. Tracking of Nonstatio-nary EEG with the Roots of ARMA Models. In IEEE Conf. EMBC-95, 1995. Dis-ponıvel em http://citeseer.ist.psu.edu/patomaki95tracking.html, aces-sado em 09/12/2004.

[33] PATOMAKI, L.; KAIPIO, J. P.; KARJALAINEN, P. A.; JUNTUNEN, M.Tracking of Nonstationary EEG with the Polynomial Root Perturbation. InIEEE Conf. EMBS-96, 1996. Disponıvel em http://citeseer.ist.psu.edu/

patomaki96tracking.html, acessado em 09/12/2004.

[34] PRESS, William H.; FLANNERY, Brian P.; TEUKOLSKY, Saul A.; VETTER-LING, William T. Numerical Recipes in C: The Art of Scientific Computing. 2nd

ed. Cambridge University Press, Cambridge, 1992.

[35] SIMON, Claudio; SHEN, Jian; SEIP, Ralf; EBBINI, Emad S. A Robust and Com-putationally Efficient Algorithm for Mean Scatterer Spacing Estimation IEEETransactions on Ultrasonics, Ferroelectrics and Frequency Control, v.44, n.4,p.882-894, 1997.

[36] SOBOL, I. Iniciacao na Matematica - O Metodo de Monte Carlo. Ed. Mir Moscou,1983.

135


[37] SODSRI, CHUKIET. Time-Varying Autoregressive Modelling for NonstationaryAcoustic Signal and Its Frequency Analysis. The Pennsylvania State University,Thesis in Acoustics, USA, 2003.

[38] SCALASSARA, Paulo Rogerio. Medidas de Caracterısticas Fısicas de TecidosBiologicos por Inspecao de Ultra-som. Monografia (Trabalho de Conclusao deCurso em Engenharia Eletrica). Orientador: Carlos Dias Maciel. UniversidadeEstadual de Londrina, Londrina, 2003.

[39] SCALASSARA, Paulo Rogerio; NEGRAO, Fernando Alves; TSUZUKI, EdgarMatsuo; MACIEL, Carlos Dias. Medida de Velocidade de Propagacao do Som eEspessura de Materiais por Ultra-som. SEMINA: Ciencias Exatas e Tecnologicas,v.23, p.9-16, Londrina, 2002.

[40] TARVAINEN, Mika. Estimation Methods for Nonstationary Biosignals. Universityof Kuopio, Department of Applied Physics, Doctoral Dissertation, Finland, 2004.

[41] THERRIEN, Charles W. Discrete Random Signals and Statistical Signal Pro-cessing. Prentice Hall Signal Processing Series. Prentice Hall, Inc., New Jersey,1992.

[42] TSUZUKI, Edgar Matsuo. Desenvolvimento de Software para Aquisicao de Da-dos de Ultra-som. Monografia (Trabalho de Conclusao de Curso em EngenhariaEletrica). Orientador: Carlos Dias Maciel. Universidade Estadual de Londrina,Londrina, 2003.

[43] VARGA, Andras. A Multishift Hessenberg Method for Pole Assignment of Single-Input Systems. IEEE Transactions on Automatic Control, v.41, n.12, p.1795-1799,1996.

[44] VARGHESE, Tomy; DONOHUE, Kevin D.; GENIS, Vlad I.; HALPERN, EthanJ. Order Selection Criteria for Detecting Mean Scatterer Spacings with the ARModel. IEEE Transactions on Ultrasonics, Ferroelectrics and Frequency Control,v.43, n.5, p. 979-984, 1996.

[45] VARGHESE, Tomy; DONOHUE, Kevin D. Estimating Mean Scatterer Spacingwith the Frequency-Smoothed Spectral Autocorrelation Function. IEEE Transac-tions on Ultrasonics, Ferroelectrics and Frequency Control, v.42, n.3, p.451-463,1995.

[46] VARGHESE, Tomy; DONOHUE, Kevin D. Mean Scatterer Spacing Estimateswith Spectral Correlation. Journal of The Acoustical Society of America, v.96,n.6, p.3504-3514, 1994.

[47] WANG, Tao; SANIIE, Jafar; JIN, Xiaomei. Analysis of Low-Order Autoregres-sive Models for Ultrasonic Grain Signal Characterization. IEEE Transactions onUltrasonics, Ferroelectrics and Frequency Control, v.38, n.2, p.116-124, 1991.

136


[48] WEAR, K. A.; WAGNER, R. F.; BRIAN, S. G. A Comparison of AutoregressiveSpectral Estimation Algorithms and Order Determination Methods in UltrasonicTissue Characterization. IEEE Transactions on Ultrasonics, Ferroelectrics andFrequency Control, v.42, n.4, p. 709-716, 1995.

[49] WEAR, Keith A.; WAGNER, Robert F.; INSANA, Michael F.; HALL, TimothyJ. Application of Autoregressive Spectral Analysis to Cepstral Estimation of MeanScatterer Spacing. IEEE Transactions on Ultrasonics, Ferroelectrics and FrequencyControl, v.40, n.1, p. 50-58, 1993.

[50] WEBSTER, J. G. Medical Instrumentation: Application and Design. 2nd ed. JohnWiley & Sons, New York, 1995.

137

Documents

An´alise de Sinais de Ultra-som usando Decomposic¸ao …la2i.com/meel/scalassara05analise.pdf · - Um outro agradecimento a Erica, bolsista junior´ do CNPq, pelo aux´ılio na