New upstream version 6.0.52 (e5c9446e) · Commits · Debian Med / praat

EEG/EEG.cpp

+56 −44

Original line number	Diff line number	Diff line
		/* EEG.cpp
		*
		* Copyright (C) 2011-2018 Paul Boersma
		* Copyright (C) 2011-2019 Paul Boersma
		*
		* This code is free software; you can redistribute it and/or modify
		* it under the terms of the GNU General Public License as published by
		@@ -97,28 +97,38 @@ autoEEG EEG_readFromBdfFile (MelderFile file) {
		try {
		autofile f = Melder_fopen (file, "rb");
		char buffer [81];
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		bool is24bit = buffer [0] == (char) 255;
		fread (buffer, 1, 80, f); buffer [80] = '\0';
		fread (buffer, 1, 80, f);
		buffer [80] = '\0';
		trace (U"Local subject identification: \"", Melder_peek8to32 (buffer), U"\"");
		fread (buffer, 1, 80, f); buffer [80] = '\0';
		fread (buffer, 1, 80, f);
		buffer [80] = '\0';
		trace (U"Local recording identification: \"", Melder_peek8to32 (buffer), U"\"");
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		trace (U"Start date of recording: \"", Melder_peek8to32 (buffer), U"\"");
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		trace (U"Start time of recording: \"", Melder_peek8to32 (buffer), U"\"");
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		integer numberOfBytesInHeaderRecord = atol (buffer);
		trace (U"Number of bytes in header record: ", numberOfBytesInHeaderRecord);
		fread (buffer, 1, 44, f); buffer [44] = '\0';
		fread (buffer, 1, 44, f);
		buffer [44] = '\0';
		trace (U"Version of data format: \"", Melder_peek8to32 (buffer), U"\"");
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		integer numberOfDataRecords = strtol (buffer, nullptr, 10);
		trace (U"Number of data records: ", numberOfDataRecords);
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		double durationOfDataRecord = atof (buffer);
		trace (U"Duration of a data record: ", durationOfDataRecord);
		fread (buffer, 1, 4, f); buffer [4] = '\0';
		fread (buffer, 1, 4, f);
		buffer [4] = '\0';
		integer numberOfChannels = atol (buffer);
		trace (U"Number of channels in data record: ", numberOfChannels);
		if (numberOfBytesInHeaderRecord != (numberOfChannels + 1) * 256)
		@@ -126,7 +136,8 @@ autoEEG EEG_readFromBdfFile (MelderFile file) {
		U") doesn't match number of channels (", numberOfChannels, U").");
		autostring32vector channelNames (numberOfChannels);
		for (integer ichannel = 1; ichannel <= numberOfChannels; ichannel ++) {
		fread (buffer, 1, 16, f); buffer [16] = '\0'; // labels of the channels
		fread (buffer, 1, 16, f);
		buffer [16] = '\0'; // labels of the channels
		/*
		* Strip all final spaces.
		*/
		@@ -143,37 +154,45 @@ autoEEG EEG_readFromBdfFile (MelderFile file) {
		bool hasLetters = str32equ (channelNames [numberOfChannels].get(), U"EDF Annotations");
		double samplingFrequency = undefined;
		for (integer channel = 1; channel <= numberOfChannels; channel ++) {
		fread (buffer, 1, 80, f); buffer [80] = '\0'; // transducer type
		fread (buffer, 1, 80, f);
		buffer [80] = '\0'; // transducer type
		}
		for (integer channel = 1; channel <= numberOfChannels; channel ++) {
		fread (buffer, 1, 8, f); buffer [8] = '\0'; // physical dimension of channels
		fread (buffer, 1, 8, f);
		buffer [8] = '\0'; // physical dimension of channels
		}
		autoVEC physicalMinimum (numberOfChannels, kTensorInitializationType::RAW);
		for (integer ichannel = 1; ichannel <= numberOfChannels; ichannel ++) {
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		physicalMinimum [ichannel] = atof (buffer);
		}
		autoVEC physicalMaximum (numberOfChannels, kTensorInitializationType::RAW);
		for (integer ichannel = 1; ichannel <= numberOfChannels; ichannel ++) {
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		physicalMaximum [ichannel] = atof (buffer);
		}
		autoVEC digitalMinimum (numberOfChannels, kTensorInitializationType::RAW);
		for (integer ichannel = 1; ichannel <= numberOfChannels; ichannel ++) {
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		digitalMinimum [ichannel] = atof (buffer);
		}
		autoVEC digitalMaximum (numberOfChannels, kTensorInitializationType::RAW);
		for (integer ichannel = 1; ichannel <= numberOfChannels; ichannel ++) {
		fread (buffer, 1, 8, f); buffer [8] = '\0';
		fread (buffer, 1, 8, f);
		buffer [8] = '\0';
		digitalMaximum [ichannel] = atof (buffer);
		}
		for (integer channel = 1; channel <= numberOfChannels; channel ++) {
		fread (buffer, 1, 80, f); buffer [80] = '\0'; // prefiltering
		fread (buffer, 1, 80, f);
		buffer [80] = '\0'; // prefiltering
		}
		integer numberOfSamplesPerDataRecord = 0;
		for (integer channel = 1; channel <= numberOfChannels; channel ++) {
		fread (buffer, 1, 8, f); buffer [8] = '\0'; // number of samples in each data record
		fread (buffer, 1, 8, f);
		buffer [8] = '\0'; // number of samples in each data record
		integer numberOfSamplesInThisDataRecord = atol (buffer);
		if (isundef (samplingFrequency)) {
		numberOfSamplesPerDataRecord = numberOfSamplesInThisDataRecord;
		@@ -185,7 +204,8 @@ autoEEG EEG_readFromBdfFile (MelderFile file) {
		U") doesn't match sampling frequency of channel 1 (", samplingFrequency, U").");
		}
		for (integer channel = 1; channel <= numberOfChannels; channel ++) {
		fread (buffer, 1, 32, f); buffer [32] = '\0'; // reserved
		fread (buffer, 1, 32, f);
		buffer [32] = '\0'; // reserved
		}
		double duration = numberOfDataRecords * durationOfDataRecord;
		autoEEG him = EEG_create (0, duration);
		@@ -282,8 +302,10 @@ autoEEG EEG_readFromBdfFile (MelderFile file) {
		time = undefined; // defensive
		}
		} else {
		thee = TextGrid_create (0, duration,
		numberOfStatusBits == 8 ? U"S1 S2 S3 S4 S5 S6 S7 S8" : U"S1 S2 S3 S4 S5 S6 S7 S8 S9 S10 S11 S12 S13 S14 S15 S16", U"");
		thee = TextGrid_create (0.0, duration,
		numberOfStatusBits == 8 ? U"S1 S2 S3 S4 S5 S6 S7 S8" : U"S1 S2 S3 S4 S5 S6 S7 S8 S9 S10 S11 S12 S13 S14 S15 S16",
		U""
		);
		for (int bit = 1; bit <= numberOfStatusBits; bit ++) {
		uint32 bitValue = 1 << (bit - 1);
		IntervalTier tier = (IntervalTier) thy tiers->at [bit];
		@@ -435,11 +457,10 @@ void EEG_filter (EEG me, double lowFrequency, double lowWidth, double highFreque
		autoSpectrum spec = Sound_to_Spectrum (channel.get(), true);
		Spectrum_passHannBand (spec.get(), lowFrequency, 0.0, lowWidth);
		Spectrum_passHannBand (spec.get(), 0.0, highFrequency, highWidth);
		if (doNotch50Hz) {
		if (doNotch50Hz)
		Spectrum_stopHannBand (spec.get(), 48.0, 52.0, 1.0);
		}
		autoSound him = Spectrum_to_Sound (spec.get());
		NUMvector_copyElements (& his z [1] [0], & my sound -> z [ichan] [0], 1, my sound -> nx);
		my sound -> z.row (ichan) <<= his z.row (1).part (1, my sound -> nx);
		}
		} catch (MelderError) {
		Melder_throw (me, U": not filtered.");
		@@ -472,15 +493,13 @@ void EEG_subtractReference (EEG me, conststring32 channelName1, conststring32 ch
		if (channelNumber1 == 0)
		Melder_throw (me, U": no channel named \"", channelName1, U"\".");
		integer channelNumber2 = EEG_getChannelNumber (me, channelName2);
		if (channelNumber2 == 0 && channelName2 [0] != '\0')
		if (channelNumber2 == 0 && channelName2 [0] != U'\0')
		Melder_throw (me, U": no channel named \"", channelName2, U"\".");
		const integer numberOfElectrodeChannels = my numberOfChannels - EEG_getNumberOfExtraSensors (me);
		for (integer isamp = 1; isamp <= my sound -> nx; isamp ++) {
		double referenceValue = channelNumber2 == 0 ? my sound -> z [channelNumber1] [isamp] :
		0.5 * (my sound -> z [channelNumber1] [isamp] + my sound -> z [channelNumber2] [isamp]);
		for (integer ichan = 1; ichan <= numberOfElectrodeChannels; ichan ++) {
		my sound -> z [ichan] [isamp] -= referenceValue;
		}
		const double referenceValue = ( channelNumber2 == 0 ? my sound -> z [channelNumber1] [isamp] :
		0.5 * (my sound -> z [channelNumber1] [isamp] + my sound -> z [channelNumber2] [isamp]) );
		my sound -> z.column (isamp).part (1, numberOfElectrodeChannels) -= referenceValue;
		}
		}

		@@ -493,12 +512,8 @@ void EEG_subtractMeanChannel (EEG me, integer fromChannel, integer toChannel) {
		Melder_throw (U"Channel range cannot run from ", fromChannel, U" to ", toChannel, U". Please reverse.");
		const integer numberOfElectrodeChannels = my numberOfChannels - EEG_getNumberOfExtraSensors (me);
		for (integer isamp = 1; isamp <= my sound -> nx; isamp ++) {
		double referenceValue = 0.0;
		for (integer ichan = fromChannel; ichan <= toChannel; ichan ++)
		referenceValue += my sound -> z [ichan] [isamp];
		referenceValue /= (toChannel - fromChannel + 1);
		for (integer ichan = 1; ichan <= numberOfElectrodeChannels; ichan ++)
		my sound -> z [ichan] [isamp] -= referenceValue;
		const double referenceValue = NUMmean (my sound -> z.column (isamp).part (fromChannel, toChannel));
		my sound -> z.column (isamp).part (1, numberOfElectrodeChannels) -= referenceValue;
		}
		}

		@@ -506,10 +521,7 @@ void EEG_setChannelToZero (EEG me, integer channelNumber) {
		try {
		if (channelNumber < 1 \|\| channelNumber > my numberOfChannels)
		Melder_throw (U"No channel ", channelNumber, U".");
		integer numberOfSamples = my sound -> nx;
		VEC channel = my sound -> z.row (channelNumber);
		for (integer isample = 1; isample <= numberOfSamples; isample ++)
		channel [isample] = 0.0;
		my sound -> z.row (channelNumber) <<= 0.0;
		} catch (MelderError) {
		Melder_throw (me, U": channel ", channelNumber, U" not set to zero.");
		}
		@@ -563,7 +575,7 @@ autoEEG EEG_extractChannel (EEG me, conststring32 channelName) {
		}
		}

		autoEEG EEG_extractChannels (EEG me, constVEC channelNumbers) {
		autoEEG EEG_extractChannels (EEG me, constVECVU const& channelNumbers) {
		try {
		integer numberOfChannels = channelNumbers.size;
		Melder_require (numberOfChannels > 0,
		@@ -590,7 +602,7 @@ static void Sound_removeChannel (Sound me, integer channelNumber) {
		Melder_require (my ny > 1,
		U"Cannot remove last remaining channel.");
		for (integer ichan = channelNumber; ichan < my ny; ichan ++)
		NUMvector_copyElements (& my z [ichan + 1] [0], & my z [ichan] [0], 1, my nx);
		my z.row (ichan) <<= my z.row (ichan + 1);
		my ymax -= 1.0;
		my ny -= 1;
		} catch (MelderError) {

EEG/EEG.h

+1 −1

Original line number	Diff line number	Diff line
		@@ -53,7 +53,7 @@ void EEG_setChannelToZero (EEG me, conststring32 channelName);
		void EEG_removeTriggers (EEG me, kMelder_string which, conststring32 criterion);
		autoEEG EEG_extractChannel (EEG me, integer channelNumber);
		autoEEG EEG_extractChannel (EEG me, conststring32 channelName);
		autoEEG EEG_extractChannels (EEG me, constVEC channelNumbers);
		autoEEG EEG_extractChannels (EEG me, constVECVU const& channelNumbers);
		void EEG_removeChannel (EEG me, integer channelNumber);
		void EEG_removeChannel (EEG me, conststring32 channelName);
		static inline autoSound EEG_extractSound (EEG me) { return Data_copy (my sound.get()); }

LPC/Cepstrogram.cpp

+2 −6

Original line number	Diff line number	Diff line
		@@ -184,22 +184,18 @@ autoPowerCepstrogram PowerCepstrogram_smooth (PowerCepstrogram me, double timeAv
		integer numberOfFrames = Melder_ifloor (timeAveragingWindow / my dx);
		if (numberOfFrames > 1) {
		autoVEC qin = newVECraw (my nx);
		autoVEC qout = newVECraw (my nx);
		for (integer iq = 1; iq <= my ny; iq ++) {
		qin.all() <<= thy z.row (iq); // ppgb: why this extra copying?
		VECsmoothByMovingAverage_preallocated (qout.get(), qin.get(), numberOfFrames);
		thy z.row (iq) <<= qout.all();
		VECsmoothByMovingAverage_preallocated (thy z.row (iq), qin.get(), numberOfFrames);
		}
		}
		// 2. average across quefrencies
		integer numberOfQuefrencyBins = Melder_ifloor (quefrencyAveragingWindow / my dy);
		if (numberOfQuefrencyBins > 1) {
		autoVEC qin = newVECraw (thy ny);
		autoVEC qout = newVECraw (thy ny);
		for (integer iframe = 1; iframe <= my nx; iframe ++) {
		qin.get() <<= thy z.column (iframe);
		VECsmoothByMovingAverage_preallocated (qout.get(), qin.get(), numberOfQuefrencyBins);
		thy z.column (iframe) <<= qout.get();
		VECsmoothByMovingAverage_preallocated (thy z.column (iframe), qin.get(), numberOfQuefrencyBins);
		}
		}
		return thee;

LPC/Cepstrum.cpp

+1 −1

Original line number	Diff line number	Diff line
		@@ -110,7 +110,7 @@ static void _Cepstrum_draw (Cepstrum me, Graphics g, double qmin, double qmax, d
		if (autoscaling)
		NUMextrema (y.get(), & minimum, & maximum);
		else
		VECclip_inplace (y.get(), minimum, maximum);
		VECclip_inplace_inline (y.get(), minimum, maximum);

		Graphics_setWindow (g, qmin, qmax, minimum, maximum);
		Graphics_function (g, y.at, 1, numberOfSelected, Matrix_columnToX (me, imin), Matrix_columnToX (me, imax));

LPC/LPC.cpp

+1 −1

Original line number	Diff line number	Diff line
		@@ -125,7 +125,7 @@ void LPC_drawGain (LPC me, Graphics g, double tmin, double tmax, double gmin, do
		void LPC_drawPoles (LPC me, Graphics g, double time, bool garnish) {
		autoPolynomial p = LPC_to_Polynomial (me, time);
		autoRoots r = Polynomial_to_Roots (p.get());
		Roots_draw (r.get(), g, -1.0, 1.0, -1.0, 1.0, U"+", 12, garnish);
		Roots_draw (r.get(), g, -1.0, 1.0, -1.0, 1.0, U"+", 12.0, garnish);
		}

		autoMatrix LPC_downto_Matrix_lpc (LPC me) {