Repository navigation
Expand file tree
/
Copy pathCLASS_codec.m
More file actions
2883 lines (2535 loc) · 145 KB
/
Copy pathCLASS_codec.m
File metadata and controls
2883 lines (2535 loc) · 145 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
%> @file CLASS_codec.cpp
% ======================================================================
%> @brief CLASS_codec encapuslates much of the coding and decoding
%> functionality required by SEV's IO routines. The methods here are all static.
%> @note Author: Hyatt Moore IV
%> @note Created 9/29/2014 - Derived from CLASS_events_container static
%> methods.
% ======================================================================
classdef CLASS_codec < handle
properties(Constant)
SECONDS_PER_EPOCH = 30;
EDF_ANNOTATIONS_CHANNEL = 1;
end
methods(Static)
% ======================================================================
%> @brief Loads staging data from hypnogram file, applies unknown
%> staging label (Default = 7) as applicable, determines sleep
%> cycle period, and extract first non wake epoch.
% ======================================================================
%> @param stages_filename is the filename of an ASCII tab-delimited file whose
%> second column contains a vector of scored sleep stages for each epoch of
%> a sleep study.
%> @param num_epochs
%> @param unknown_stage_label (Optional) Integer number to use for
%> unclassified/unknown hypnogram stages. Default is 7.
%> @retval STAGES A struct with the following fields
%> - @c line = the second column of stages_filename - the scored sleep stages
%> - @c count = the number of stages for each one
%> - @c cycle - the nrem/rem cycle
%> - @c firstNonWake - index of first non-Wake(0) and non-unknown(7) staged epoch
%> - @c standard_epoch_sec 30 seconds
%> - @c filename Name of the staging filename
%> - @c study_duration_in_seconds How long the study is measured in
%> seconds based on the number of epochs entered and the
%> standard_epoch_sec duration (i.e. 30 seconds)
%> - @c standard_epoch_sec standard epoch duration in sec (30 s)
%> - @c standard_epoch_min minutes per epoch (e.g. 0.5 m)
%> @note Author: Hyatt Moore IV
%> Written: 9.26.2012
%> modified before 12.3.2012 to include scoreSleepCycles(.);
%> modified 1/16/2013 - added .firstNonWake
%> modified 2/2/2013 - added .standard_epoch_sec = 30
%> .study_duration_in_seconds
%> modified 5.1.2013 - added .filename = stages_filename;
%> modified 12.17.2019 - adding stage summary metrics (e.g. sleep
%> latency metric)
function STAGES = loadSTAGES(stages_filename,num_epochs,unknown_stage_label)
if(nargin<3)
default_unknown_stage = 7;
if(nargin<2)
num_epochs = -1;
if(nargin<1)
stages_filename = [];
end
end
else
default_unknown_stage = unknown_stage_label;
end
% Make a default line, regardless of what comes out of here.
stages = repmat(default_unknown_stage,num_epochs,1);
if(~isempty(stages_filename))
[~,~,ext] = fileparts(stages_filename);
%load stages information if the file exists and we know its
%extension.
if(exist(stages_filename,'file') && (strcmpi(ext,'.sta')||strcmpi(ext,'.evts')))
if(strcmpi(ext,'.sta') || strcmpi(ext,'.evts'))
if(strcmpi(ext,'.sta'))
stages = load(stages_filename,'-ASCII'); %for ASCII file type loading
else
[~,stages] = CLASS_codec.parseSSCevtsFile(stages_filename,default_unknown_stage);
end
% elseif(strcmpi(ext,'.evts'))
% [~,stageVec] = CLASS_codec.parseSSCevtsFile(stages_filename,default_unknown_stage);
% STAGES.line = stageVec;
% STAGES.epochs = 1:numel(STAGES.line);
% end
end
else
mfile = strcat(mfilename('fullpath'),'.m');
fprintf('failed in %s\n\tFilename argument for loadSTAGES could not be found.\n',mfile);
% throw(MException('SEV:ARGERR','Filename argument for loadSTAGES could not be found'));
end
else
mfile = strcat(mfilename('fullpath'),'.m');
fprintf('failed in %s\n\tMissing or empty filename argument for loadSTAGES\n',mfile);
% throw(MException('SEV:ARGERR','Missing or empty filename argument for loadSTAGES'));
end
STAGES = CLASS_codec.stages2STAGES(stages,num_epochs, default_unknown_stage);
STAGES.filename = stages_filename;
end
% Returns axes handle (ax)
function ax = plotSTAGES(STAGES, study_start_time)
%
% W (w) 0
% R (r) 5
% N1 (p) 1
% N2 (v) 2
% N3 (g) 3
% Lights out - (next hour skipped) - then each hour (HH AM/PM)
if nargin<2 || isempty(study_start_time)
study_start_time = [0 0 0];
end
hypno = STAGES.line(STAGES.tib_first_epoch:STAGES.tib_final_epoch);
% lights out
staging_start_time = datenum([0 0 0, study_start_time]+[0 0 0 0 STAGES.tib_first_epoch*0.5 0]);
epoch_time = datenum([0 0 0 0 0 30]);
% lights on
staging_end_time = datenum([0 0 0, study_start_time]+[0 0 0 0 STAGES.tib_final_epoch*0.5 0]);
% end_time = start_time+(numel(hypno)-1)*epoch_time;
x = staging_start_time:epoch_time:staging_end_time;
% go every hour after an hour, rounded up, from the start time.
x_ticks = [staging_start_time, ceil(mod(staging_start_time, 1)*24+1)/24:1/24:staging_end_time];
hypno_offset = [0 10 'k'+0
5 8 'r'+0
1 6 'm'+0
2 4 'b'+0
3 2 'g'+0
4 2 'g'+0
7 0 'w'+0] ;
%hypno_line = nan(size(hypno));
% Place the lines at the specified offsets.
set(gcf,'position', [10 1000 750 128]);
ax = axes('box','on', 'nextplot','add','tickdir','out',...
'ylim',[0, 12], 'xlim',[staging_start_time, staging_end_time],...
'xtickmode','manual','xtick',x_ticks,...
'plotboxaspectratiomode','manual','plotboxaspectratio',[30 12 1],...
'dataaspectratiomode','manual','dataaspectratio',[3 1200 1],...
'yticklabelmode','manual','ytick', flipud(hypno_offset(1:end-2,2)),...
'yticklabel',fliplr({'W','R','N1','N2','N3'}));
width = 1;
for s=1:size(hypno_offset,1)
hypno_line = nan(size(hypno));
hypno_line(hypno==hypno_offset(s,1)) = hypno_offset(s,2);
bar(ax, x, hypno_line, width, char(hypno_offset(s,3)));
end
% datetick(ax,'x','HH PM','keeplimits');
datetick(ax,'x','HH PM','keepticks','keeplimits');
xticklabels = string(get(ax,'xticklabel'));
xticklabels{1} = datestr(staging_start_time,'HH:MM:SS PM');
set(ax,'xticklabel',xticklabels);
end
% stages is either an Mx1 vector or an Mx2 or Mx3 matrix
% If Mx1 then elements are stages in epoch order.
% If Mx2 then column 1 comprises the scored epoch number and column 2
% comprises the sleep stage scored for the corresponding epoch
function STAGES = stages2STAGES(stages, num_epochs, default_unknown_stage)
narginchk(1,3);
if(nargin<3)
default_unknown_stage = 7;
end
STAGES.standard_epoch_sec = CLASS_codec.SECONDS_PER_EPOCH; %30 second epochs
STAGES.standard_epoch_min = CLASS_codec.SECONDS_PER_EPOCH/60; %30 second epochs
if(size(stages,2)>2)
stages = stages(:, 1:2);
end
if(size(stages,2)~=2)
stage_vec = stages;
epochs = 1:numel(stage_vec);
stages = [epochs(:), stage_vec];
end
% stages is Mx2 matrix with column 1 being the
% epochs and column 2 being the stage
if(nargin<2 || num_epochs<0)
num_epochs = size(stages,1);
end
if(nargin>1 && ~isempty(num_epochs) && floor(num_epochs)>0)
if(num_epochs~=size(stages,1))
STAGES.epochs = stages(:,1);
STAGES.line = repmat(default_unknown_stage,max([num_epochs;size(stages,1);STAGES.epochs(:)]),1);
STAGES.line(STAGES.epochs) = stages(:,2);
else
%this cuts things off at the end, where we assume the
%disconnect between num_epochs expected and num epochs found
%has occurred. However, logically, there is no guarantee that
%the disconnect did not occur anywhere else (e.g. at the
%beginning, or sporadically throughout)
STAGES.line = stages(1:floor(num_epochs),2);
end
else
STAGES.line = stages(:,2); %grab the sleep stages
end
STAGES.filename = [];
STAGES.firstNonWake = [];
STAGES.summary_field_names = CLASS_codec.getSTAGESSummaryFieldNames();
STAGES.count = zeros(8,1);
STAGES.bouts = zeros(8,1);
if(nargin<2)
num_epochs = numel(STAGES.line);
end
%calculate number of epochs in each stage
for k = 0:numel(STAGES.count)-1
STAGES.count(k+1) = sum(STAGES.line==k);
bout_diff = diff(STAGES.line==k);
STAGES.bouts(k+1) = max(sum(bout_diff==1), sum(bout_diff==-1));
end
firstNonWake = 1;
while( firstNonWake<=numel(STAGES.line) && (STAGES.line(firstNonWake)==7||STAGES.line(firstNonWake)==0))
firstNonWake = firstNonWake+1;
end
if firstNonWake> numel(STAGES.line)
firstNonWake = [];
end
STAGES.firstN1_epoch = find(STAGES.line==1,1);
STAGES.firstN2_epoch = find(STAGES.line==2,1);
STAGES.firstSWS_epoch = find(STAGES.line==3|STAGES.line==4,1);
[STAGES.tib_first_epoch, STAGES.tib_final_epoch] = CLASS_codec.getTIBEpochs(STAGES.line, default_unknown_stage);
lastNonWake = STAGES.tib_final_epoch;
if ~isempty(lastNonWake)
while( lastNonWake>=1 && (STAGES.line(lastNonWake)==7||STAGES.line(lastNonWake)==0))
lastNonWake = lastNonWake-1;
end
if lastNonWake<1
lastNonWake = [];
end
end
STAGES.firstNonWake = firstNonWake;
STAGES.lastNonWake = lastNonWake;
if(num_epochs~=numel(STAGES.line))
fprintf(1,'Hypnogram contains %u stages, but shows it should have %u\n', numel(STAGES.line),num_epochs);
end
% tib
wk_tib_bout_diff = diff(STAGES.line(STAGES.tib_first_epoch:STAGES.tib_final_epoch)==0);
STAGES.bouts_wake_during_tib = max(sum(wk_tib_bout_diff==1), sum(wk_tib_bout_diff==-1));
% sleep period time
wk_spt_bout_diff = diff(STAGES.line(STAGES.firstNonWake:STAGES.lastNonWake)==0);
STAGES.bouts_wake_during_spt = max(sum(wk_spt_bout_diff==1), sum(wk_spt_bout_diff==-1));
%this may be unnecessary when the user does not care about sleep cycles.
% STAGES.cycles = scoreSleepCycles(STAGES.line);
STAGES.cycles = scoreSleepCycles_ver_REMweight(STAGES.line);
STAGES.study_duration_in_seconds = STAGES.standard_epoch_sec*numel(STAGES.line);
STAGES.study_duration_minutes = numel(STAGES.line) * STAGES.standard_epoch_min;
STAGES.total_wake_sleep_minutes = sum(STAGES.line~=7) * STAGES.standard_epoch_min;
STAGES.sleep_latency_minutes = CLASS_codec.getSleepLatency(STAGES);
STAGES.rem_latency_minutes = CLASS_codec.getREMLatency(STAGES);
STAGES.n1_latency_minutes = nan;
STAGES.n2_latency_minutes = nan;
STAGES.sws_latency_minutes = nan;
if ~isempty(STAGES.firstN1_epoch)
STAGES.n1_latency_minutes = (STAGES.firstN1_epoch-STAGES.tib_first_epoch)*STAGES.standard_epoch_min;
end
if ~isempty(STAGES.firstN2_epoch)
STAGES.n2_latency_minutes = (STAGES.firstN2_epoch-STAGES.tib_first_epoch)*STAGES.standard_epoch_min;
end
if ~isempty(STAGES.firstSWS_epoch)
STAGES.sws_latency_minutes = (STAGES.firstSWS_epoch-STAGES.tib_first_epoch)*STAGES.standard_epoch_min;
end
STAGES.waso_minutes = sum(STAGES.line(STAGES.firstNonWake:end)==0) * STAGES.standard_epoch_min;
STAGES.total_wake_min = sum(STAGES.line==0) * STAGES.standard_epoch_min;
STAGES.total_sleep_min = sum(STAGES.line~=0 & STAGES.line~=7) * STAGES.standard_epoch_min;
STAGES.total_s1_min = sum(STAGES.line==1) * STAGES.standard_epoch_min;
STAGES.total_s2_min = sum(STAGES.line==2) * STAGES.standard_epoch_min;
STAGES.total_sws_min = sum(STAGES.line==3|STAGES.line==4) * STAGES.standard_epoch_min;
STAGES.total_rem_min = sum(STAGES.line==5) * STAGES.standard_epoch_min;
% spt - sleep period time
STAGES.total_sleep_period_min = (STAGES.lastNonWake-STAGES.firstNonWake+1) * STAGES.standard_epoch_min;
STAGES.total_wake_during_tib_min = sum(STAGES.line(STAGES.tib_first_epoch: STAGES.tib_final_epoch)==0) * STAGES.standard_epoch_min;
STAGES.total_wake_during_sleep_period_min = sum(STAGES.line(STAGES.firstNonWake: STAGES.lastNonWake)==0) * STAGES.standard_epoch_min;
num_epochs_in_bed = STAGES.tib_final_epoch-STAGES.tib_first_epoch+1;
if isempty(num_epochs_in_bed)
num_epochs_in_bed = 0;
end
STAGES.tib_min = max(num_epochs_in_bed *STAGES.standard_epoch_min,0); % ensure sure we don't get negative values with max(...,0)
if STAGES.tib_min > 0
STAGES.sleep_efficiency = STAGES.total_sleep_min/STAGES.tib_min*100;
else
STAGES.sleep_efficiency = 0;
end
if STAGES.total_sleep_min == 0
STAGES.percent_s1 = 0;
STAGES.percent_s2 = 0;
STAGES.percent_sws = 0;
STAGES.percent_rem = 0;
else
STAGES.percent_s1 = STAGES.total_s1_min / STAGES.total_sleep_min*100;
STAGES.percent_s2 = STAGES.total_s2_min / STAGES.total_sleep_min*100;
STAGES.percent_sws = STAGES.total_sws_min / STAGES.total_sleep_min*100;
STAGES.percent_rem = STAGES.total_rem_min / STAGES.total_sleep_min*100;
end
% Duration in minutes of epochs scored as 7 ('unknown')
STAGES.unscored_tib_min = STAGES.tib_min - STAGES.total_wake_sleep_minutes;
if (STAGES.unscored_tib_min < 0)
%A0120_6
fprintf(1,'Warning\n');
end
end
function [first_epoch, final_epoch] = getTIBEpochs(stages_vec, not_in_bed_value)
if nargin<2
not_in_bed_value = 7;
end
first_epoch = 1;
final_epoch = numel(stages_vec);
mask_indices = stages_vec==not_in_bed_value;
if(mask_indices(1))
first_epoch = find(~mask_indices, 1, 'first');
end
% book end - does it end masked (i.e. finish with 7's)
if(mask_indices(end))
final_epoch = find(~mask_indices,1, 'last');
end
end
function summary_field_names = getSTAGESSummaryFieldNames()
summary_field_names = { ...
'study_duration_minutes', 'tib_min', ...
'total_wake_sleep_minutes', 'unscored_tib_min', ...
'sleep_latency_minutes', 'rem_latency_minutes', 'waso_minutes', ...
'total_wake_min', 'total_sleep_min', 'total_s1_min', ...
'total_s2_min', 'total_sws_min', 'total_rem_min' ...
};
end
function sleepl_minutes_elapsed = getSleepLatency(STAGES)
sleepl_minutes_elapsed = (STAGES.firstNonWake-STAGES.tib_first_epoch)*STAGES.standard_epoch_min;
end
% The presence of a SOREMP upon first falling asleep was first reported by
% Vogel in 1960 and Rechaffen and Dement in approximately half of untreated
% patients with narcolepsy. In a large study of 648 Type 1 patients and
% 1369 controls we found that a REM latency (REML) ? 15 minutes was highly
% specific in controls (98.0% but only had a specificity of 49.7%. A
% cutoff for REML of <= 15 minutes was used as it is the value also used in
% MSLT naps, although in many models REML ranging <= 20 minutes were more
% advantageous, allowing sensitivity to be raised to 60% with limited loss
% of specificity and are likely to be used in a combined model.
function [reml_minutes_elapsed, test_result] = getREMLatency(STAGES, threshold_optional)
INIT_SOREMP_THRESHOLD_MIN = 15;
REM_STAGE = 5;
if(nargin>1 && ~isempty(threshold_optional))
threshold_min = threshold_optional;
else
threshold_min = INIT_SOREMP_THRESHOLD_MIN;
end
firstREMIndex = find(STAGES.line==REM_STAGE,1);
reml_minutes_elapsed = (firstREMIndex - STAGES.firstNonWake)*STAGES.standard_epoch_min;
test_result = reml_minutes_elapsed <= threshold_min;
end
% ======================================================================
%> @brief Retrieves SEV compatible hypnogram filename (with path)
%> based on the input edf filename provided (With path).
% ======================================================================
%> @param edf_fullfilename Filename of the EDF file to find matching hypnogram file for.
%> @note The hypnogram file is expected to be in the same directory as edf_fullfilename
%> and have a file extension of '.STA' or '.evts'
%> @retval stages_filename Filename of the hypnogram with path.
%> .STA is checked first, then .evts extension is checked if .STA is
%> not found. If neither staging file type is found (.STA or .EVTS)
%> then stages_filename is returned as empty (i.e. [])
%> @retval edf_name edf filename sans pathname.
% ======================================================================
function [stages_filename, edf_name] = getStagesFilenameFromEDF(edf_fullfilename)
[edf_path,edf_name,edf_ext] = fileparts(edf_fullfilename);
stages_filename = fullfile(edf_path,strcat(edf_name,'.STA'));
if(~exist(stages_filename,'file'))
stages_filename = fullfile(edf_path,strcat(edf_name,'.EVTS'));
if(~exist(stages_filename,'file'))
stages_filename = [];
end
end
edf_name = strcat(edf_name,edf_ext);
end
% ======================================================================
%> @brief Retrieves SEV compatible events filename (with path)
%> based on the input edf filename provided (With path).
% ======================================================================
%> @param edf_fullfilename Filename of the EDF file to find matching events file for.
%> @note The events file is expected to be in the same directory as edf_fullfilename
%> and have a file extension of '.EVTS' or '.SCO'
%> @param events_filename Filename of the hypnogram with path.
%> .EVTS is checked first, then .SCO extension is checked if .STA is
%> not found. If neither staging file type is found (.STA or .EVTS)
%> then stages_filename is returned as empty (i.e. [])
%> @retval event_Filename
%> @retval edf_name edf filename sans pathname.
% ======================================================================
function [event_Filename, edf_name] = getEventsFilenameFromEDF(edf_fullfilename)
[edf_path,edf_name,edf_ext] = fileparts(edf_fullfilename);
event_Filename = fullfile(edf_path,strcat(edf_name,'.EVTS'));
if(~exist(event_Filename,'file'))
event_Filename = fullfile(edf_path,strcat(edf_name,'.SCO'));
if(~exist(event_Filename,'file'))
event_Filename = [];
end
end
edf_name = strcat(edf_name,edf_ext);
end
% =================================================================
%> @brief Retrives a function call for files in directories with a '+'
%> prefix.
%------------------------------------------------------------------%
%> @param methodName Method's name (sans path and .m)
%> @param packageName Package name (string). Supported values include
%> - @c export
%> - @c detection
%> - @c filter
%> @retval packageMethodName - Package.method name to call the method outside of its
%> directory via its package.
%> @retval fullFilename - Name of the .m file for the method with
%> pathname. Empty if no match is found.
% =================================================================
function [packageMethodName, fullFilename] = getPackageMethodName(methodName,packageName)
candidateCategories = {'export','detection','filter'};
rootPath = fileparts(mfilename('fullpath'));
if(isempty(intersect(packageName, candidateCategories)))
packageMethodName = [];
fullFilename = [];
else
fullFilename = fullfile(rootPath,strcat('+',packageName),strcat(methodName,'.m'));
if(exist(fullFilename,'file'))
packageMethodName = strcat(packageName,'.',methodName);
else
packageMethodName = [];
fullFilename = [];
end
end
end
% =================================================================
%> @brief Retrives the information file (.inf) for the package of interest.
%> Information files describe methods used in a particular toolbox.
%------------------------------------------------------------------%
%> @param packageName Package name (string). Supported values include
%> - @c export
%> - @c detection
%> - @c filter
%> @retval methodInformationFilename Full filename of the information
%> file used to describe the toolbox associated with methodCategory.
%> Information files describe methods used in a particular toolbox.
%> Empty if no match is found.
% =================================================================
function methodInformationFilename = getMethodInformationFilename(packageName)
candidateCategories = {'export','detection','filter'};
rootPath = fileparts(mfilename('fullpath'));
if(isempty(intersect(packageName, candidateCategories)))
methodInformationFilename = [];
else
methodInformationFilename = fullfile(rootPath,strcat('+',packageName),strcat(packageName,'.inf'));
if(~exist(methodInformationFilename,'file'))
methodInformationFilename = [];
end
end
end
% =================================================================
%> @brief Retrives a function call for files in directories with a '+'
%> prefix.
%------------------------------------------------------------------%
%> @param methodName Method's name (sans path and .m)
%> @param packageName Package name that method is part of (string). Supported values include
%> - @c export
%> - @c detection
%> - @c filter
%> @retval parameterStruct
%> @retval packageMethodName - method name with package prefix (e.g. export.someMethod) to call the method in
%> the package a directory outside of its category.
%> Empty if no match is found.
% =================================================================
function [parameterStruct, packageMethodName] = getMethodParameters(methodName,packageName)
candidateCategories = {'export','detection','filter'};
rootPath = fileparts(mfilename('fullpath'));
if(isempty(intersect(packageName, candidateCategories)))
packageMethodName = [];
parameterStruct = [];
else
fullFilename = fullfile(rootPath,strcat('+',packageName),strcat(methodName,'.m'));
if(exist(fullFilename,'file'))
packageMethodName = strcat(packageName,'.',methodName);
try
parameterStruct = feval(packageMethodName);
catch me
showME(me);
fprintf('An error occurred and was caught. The variable parameterStruct will be set to []\n');
parameterStruct = [];
end
else
packageMethodName = [];
parameterStruct = [];
end
end
end
%------------------------------------------------------------------%
%> @brief Parses an export package information file (.inf) and returns
%> each rows values as a struct entry.
%------------------------------------------------------------------%
%> @param exprtInfFullFilename Full filename (path and name) of the export information
%> file to parse.
%> @retval exportMethodsStruct Struct with the following fields:
%> - @c mfilename Nx1 cell of filenames of the export method
%> - @c description Nx1 cell of descriptions for each export method.
%> - @c settingsEditor Nx1 cell of the settings editor to use for each method.
%> @note One entry per non-header row parsed of .inf file
%------------------------------------------------------------------%
function exportMethodsStruct = parseExportInfFile(exportInfFullFilename)
exportMethodsStruct = [];
if(exist(exportInfFullFilename,'file'))
fid = fopen(exportInfFullFilename,'r');
scanCell = textscan(fid,'%s %s %s','headerlines',1,'delimiter',',');
fclose(fid);
exportMethodsStruct.mfilename = scanCell{1};
exportMethodsStruct.description = scanCell{2};
exportMethodsStruct.settingsEditor = scanCell{3};
end
end
function hdr = loadHDR(edfFilename)
hdr = loadHDR(edfFilename);
end
% ======================================================================
%> @brief loads/parses the .SCO file associated with an EDF.
% ======================================================================
%> @param filename Full .SCO filename (i.e. with path) to load.
%> @param dest_samplerate Sample rate to use for SCO fields (return value)
%> @param sco_samplerate The SCO file's internal samplerate used internal =
%> this the sample rate used when generating the SCO file initially.
%> Typically the samplerate is twice the highest sampling rate used for the .EDF
%> channels as shown in the .EDF header.
%> @retval SCO is a struct with the fields
%> @li .epoch - the epoch that the scored event occured in
%> @li .start_stop_matrix - the sample point that the events begin and end on
%> @li .label - the string label used to describe the event
%> @li .duration_seconds = the duration of the event in seconds (default is 1
%> second)
%> @li .start_time - events start time as a string in hour:minute:second format
%> @param
%> @note [ul,ui,uj] = unique(SCO.label); lm_epochs = SCO.epoch(3==uj,:); lm_evts = SCO.start_stop_matrix(3==uj,:);
%> @note obtains the unique epochs for the third unique event
%> @note History
%> Author: Hyatt Moore IV (< June, 2013)
%>
%> @note modified: 4/30/2014
%> 1. Added additional input arguments, dest_samplerate and sco_samplerate
%> to help with conversion. Defaults are dest_samplerate = 100 and
%> sco_samplerate = 200
%> 2. Previously the duration_sec column was only divided by 100Hz and not
%> 200Hz, so durations could be 2x's as long if the sampling rate was not
%> 100 Hz.
%> @note modified 12/15/12
%> 3. last argument to textscan has changed to a regular expression since new
%> files had extra columns at the end which are causing problems with the
%> imports.
%> @note modified: 2/6/12 -
%> 1. Some .sco files had empty lines at the top, and we need to cruise
%> through these until something hits - so I do a getl until the scan
%> works, otherwise it will crash
%> 2. Some sco files do not convert correctly due to the sample rate being
%> done at 128. The EDF is okay, but the SCO file needs to be converted
%> correctly as well. The mode sample rate from the EDF hdr is used for
%> converting
function SCO = parseSCOfile(filename,dest_samplerate, sco_samplerate)
%> filename = '/Users/hyatt4/Documents/Sleep Project/Data/Spindle_7Jun11/A0210_3 170424.SCO';
%> filename = '/Users/hyatt4/Documents/Sleep Project/Data/Spindle_7Jun11/A0097_4 174733.SCO';
% Note: EDF_filename = [filename(1:end-3),'EDF'];
% if(exist(EDF_filename,'file'))
% HDR = loadHDR(EDF_filename);
% samplerate = max(HDR.samplerate);
% conversion_factor=100/samplerate;
% else
% samplerate = 100;
% conversion_factor = 1;
% end
%
%
%
if(nargin<3 || isempty(sco_samplerate)|| sco_samplerate<0)
sco_samplerate= 200;
end
if(nargin<2 || isempty(dest_samplerate)|| dest_samplerate<0)
dest_samplerate= 100;
end
if(exist(filename,'file'))
fid = fopen(filename,'r');
% x = textscan(fid,'%f %f %f %f %s %f %s %f %*f','delimiter','\t');
% x = textscan(fid,'%f %f %f %f %s %f %s %f %*[.]','delimiter','\t');
x = textscan(fid,'%f %f %f %f %s %f %s %f %*[^\n]','delimiter','\t');
if(isempty(x{1}))
file_size_bytes = fseek(fid,0,'eof');
if(file_size_bytes==0)
disp(['File size of ',filename,' is 0. Going to crash now.']);
end
while(isempty(x{1}) && ftell(fid)<file_size_bytes)
fgetl(fid);
x = textscan(fid,'%f %f %f %f %s %f %s %f %*f','delimiter','\t');
disp(filename);
end
end
%remove potential problem of empty first line
try
if(isnan(x{1}(1)))
for k=1:numel(x)
x{k} = x{k}(2:end);
end
end
catch me
disp(me);
end
SCO.epoch = x{1};
%handle the negative case, which pops up for the more recent sco files
%which were not converted correctly the first time...
% - adjusted on 4/30/2014 reference M.Stubbs email sent August 5, 2011
neg = x{2}<0;
if(any(neg))
x{2}(neg) = abs(floor(x{2}(neg)/100));
x{3}(neg) = abs(floor(x{3}(neg)/100));
end
x{3}(x{3}==0)=300; %make the default be 1.5 second duration.
x{3}(isnan(x{3}))=300; %make the default be 1.5 second duration - typical case for Leg Movement to not have a value listed...
% EDF_filename = [filename(1:end-3),'EDF'];
% if(exist(EDF_filename,'file'))
% HDR = loadHDR(EDF_filename);
% samplerate = max(HDR.samplerate);
% conversion_factor=100/samplerate;
% else
% samplerate = 100;
% conversion_factor = 1;
% end
conversion_factor = dest_samplerate/sco_samplerate;
SCO.start_stop_matrix = floor([x{2},x{2}+x{3}]*conversion_factor); %remove any problems with the 0.5 indexing that can occur here
SCO.duration_seconds = x{3}/sco_samplerate;
SCO.start_time = x{7};
SCO.label = deblank(x{5});
fclose(fid);
else
SCO = [];
disp('filename does not exist');
end
end
%> @brief Export the current hypnogram to a .STA text file.
%> @note Hypnograms are scored at the default interval set in the
%> SEV parameters (e.g. 30 seconds)
%> @param hypnogramVec A vector of consecutively scored sleep
%> stages with numeric values interpreted as
%> - c 0 Wake
%> - c 1 Stage 1
%> - c 2 Stage 2
%> - c 3 Stage 3
%> - c 4 Stage 4
%> - c 5 REM sleep
%> - c 7 Uknown
%> @param staFilename Optional filename to save hypnogram to. If
%> not included, the user is prompted for a filename.
%> @note The .STA file format is two column, tab-delimited, ascii
%> with the first column being a consecutive list of epochs and the
%> second column containing the sleep scores that correspond to
%> the first column epochs.
function saveHypnogram2STA(hypnogramVec, staFilename)
epoch = (1:numel(hypnogramVec))';
hypnogramVec = hypnogramVec(:);
y = [epoch,hypnogramVec];
if(nargin<2 ||isempty(staFilename))
staFilename = [];
end
fid = fopen(staFilename, 'w');
if fid>2
fprintf(fid, '%d\t%d\n', y');
fclose(fid);
else
%double check that we have .STA extension
save(staFilename,'y','-ascii');
end
end
% =================================================================
%> @brief Returns the patid and study number related to the
%> patstudy key provided. This requires patstudy input string to
%> follow a particular template (see source code).
%> @param patstudy is the filename of the .edf, less the extention
%> @retval PatID
%> @retval StudyNum
% =================================================================
function [PatID,StudyNum] = getDB_PatientIdentifiers(patstudy)
%patstudy is the filename of the .edf, less the extention
% if(numel(patstudy)>=7)
% pat = '(\w{5})_(\d+)'; %WSC format
% else
% pat = '([A-Z]+)(\d+)'; %PTSD format
% end
% pat = '(SSC_\d[4]'; %SSC format
% x=regexp(patstudy,pat,'tokens');
patstudy = strrep(patstudy,'SSC_','');
x=regexp(patstudy,'(\w+)_(\d+)||(\d+)_(\d+)||([^\d]+)(\d+)','tokens');
x = x{1};
PatID = x{1};
StudyNum = x{2};
end
%> @brief stages here refers to Stanford's STAGES cohort, a multisite
%> collection of PSGs.
function [evt_Struct,samplerate_out] = parseSTAGESEventFile(eventsFile, edfFile, desired_samplerate)
samplerate_out = desired_samplerate;
if(~exist(eventsFile,'file'))
error('File does not exist or was not given: %s', eventsFile);
end
if nargin<3 || isempty(desired_samplerate)
desired_samplerate = 100;
end
if nargin<2
edfFile = [];
end
fid = fopen(eventsFile, 'r');
if fid < 1
error('Unable to open file: %s', eventsFile);
else
% Start Time,Duration (seconds),Event
% 19:52:36, 0.000, MChg
% 19:52:36, 0.000, MChg
% 19:52:36, 0.000, MChg
a = textscan(fid, '%{hh:mm:ss}T %f %s', 'headerlines', 1, 'delimiter',',');
fclose(fid);
start_time = a{1};
duration_sec = a{2};
event_label = a{3};
custom_idx = startsWith(event_label, 'custom','IgnoreCase', true);
start_time(custom_idx) = [];
duration_sec(custom_idx) = [];
event_label(custom_idx) = [];
if ~isempty(edfFile) && exist(edfFile, 'file')
HDR = loadHDR(edfFile);
else
HDR = [];
end
evt_Struct = CLASS_codec.makeEventStruct();
evt_Struct.HDR = HDR;
evt_Struct.type = eventType;
evt_Struct.start_stop_matrix = start_stop_matrix;
evt_Struct.start_sec = start_sec;
evt_Struct.stop_sec = stop_sec;
evt_Struct.dur_sec = dur_sec;
evt_Struct.epoch = epoch;
evt_Struct.stage = stage;
evt_Struct.evtID = evtID;
if(~isempty(description))
evt_Struct.description = description;
end
if(~isempty(remainder))
evt_Struct.unknown = remainder;
end
end
end
% =================================================================
%> @brief Parse events from an Embla formatted events file (.evt/.nvt)
%> The function parses Embla files
%> @li stage
%> @li plm
%> @li desat
%> @li biocals
%> @li user
%> @li numeric
%> @li baddata
%> @param evtFilename Filename (string) of the Embla formatted event file
%> @param embla_samplerate Sampling rate of the Embla event file.
%> @param desired_samplerate (optional) sample rate to convert
%> Embla events to. This is helpful when displaying a samplerate
%> different than recorded in the .evt file. If desired_samplerate
%> is not provided, then embla_samplerate is used.
%> filename in .SCO format
%> @retval embla_evt_Struct
%> @retval embla_samplerate_out
% =================================================================
function [embla_evt_Struct, embla_samplerate_out] = parseEmblaEvent(evtFilename, embla_samplerate, desired_samplerate)
%embla_samplerate_out may change if there is a difference found in
%the stage .evt file processing as determined by adjusting for
%a 30 second epoch.
if(nargin<2)
embla_samplerate = [];
end
embla_samplerate_out = embla_samplerate;
if(~exist(evtFilename,'file'))
embla_evt_Struct = [];
disp([evtFilename,' not handled']);
else
if(nargin < 3 || isempty(desired_samplerate))
desired_samplerate = 0;
end
[~,name,~] = fileparts(evtFilename);
% Get the last part then.
if contains(name, '.')
%name = extractAfter(name,'.');
name = strsplit(name,'.');
name = name{end};
end
fid = fopen(evtFilename,'r','l'); % little endian format
HDR = CLASS_codec.parseEmblaHDR(fid);
start_sec = [];
stop_sec = [];
dur_sec = [];
epoch = [];
stage = [];
start_stop_matrix = [];
description = [];
evtID = [];
remainder = [];
eventType = name;
bytes_per_uint16 = 2;
if(HDR.num_records>0 && strncmpi(deblank(HDR.label),'event',5))
fseek(fid,0,'eof');
file_size = ftell(fid);
fseek(fid,32,'bof');
bytes_remaining = file_size-ftell(fid);
bytes_per_record = bytes_remaining/HDR.num_records;
start_sample = zeros(HDR.num_records,1);
stop_sample = start_sample;
evtID = zeros(HDR.num_records,2);
if strcmpi(eventType,'arousal')
arousType = stop_sample;
arousSecondType = stop_sample;
description = cell(HDR.num_records,1);
if(false)
fseek(fid,32,'bof');
a=fread(fid,[bytes_per_record,inf],'uint8')';
a
fseek(fid,32,'bof');
end
if(bytes_per_record==52)
unknownID = zeros(HDR.num_records,3);
remainder = zeros(HDR.num_records,4,'uint8');
orderPlaced = start_sample;
for r=1:HDR.num_records
start_sample(r) = fread(fid,1,'uint32');%[0-3
stop_sample(r) = fread(fid,1,'uint32');%[4-7
evtID(r,:) = fread(fid,2,'uint8'); % [8- [6,1]
unknownID(r,:) = fread(fid,3,'uint16');%[10-15
orderPlaced(r) = fread(fid,1,'uint32');%[16-
fseek(fid,28,'cof');% [20
remainder(r,:) = fread(fid,4,'uint8');%[440
end
elseif(bytes_per_record==16)
unknownID = zeros(HDR.num_records,2);
remainder = zeros(HDR.num_records,2,'uint8');
for r=1:HDR.num_records
start_sample(r) = fread(fid,1,'uint32');
stop_sample(r) = fread(fid,1,'uint32');
evtID(r,:) = fread(fid,2,'uint8'); % [6,1]
unknownID(r,:) = fread(fid,2,'uint16');
remainder(r,:) = fread(fid,2,'uint8');
end
end
elseif(strcmpi(eventType,'baddata'))
if(false)
fseek(fid,32,'bof');
a=fread(fid,[bytes_per_record,inf],'uint8')';
fseek(fid,32,'bof');
a(:,5:end)
end
start_sample = zeros(HDR.num_records,1);
stop_sample = start_sample;
intro_size = 10;
remainder = zeros(HDR.num_records,bytes_per_record-intro_size,'uint8');
for r=1:HDR.num_records
start_sample(r) = fread(fid,1,'uint32');
stop_sample(r) = fread(fid,1,'uint32');
evtID(r,:) = fread(fid,2,'uint8'); % [5,1]
remainder(r,:) = fread(fid,bytes_per_record-intro_size,'uint8');
end
elseif(strcmpi(eventType,'biocals'))
if(false)
fseek(fid,32,'bof');
a=fread(fid,[bytes_per_record,inf],'uint8')';
a
fseek(fid,32,'bof');
end
% first line:
% [1][2] [3-4]...[23-24] [25-28] [29-32] || [33-36] [37-40] [41-42]
% [1 0] [uint16=>char] uint32 uint32 || uint32
% Title Text checksum # entries || elapsed sample start [13 1 0 0] - biocals
% [1 stage# 0 0] - stage... [34 0]
%
% Elapsed Time Format:
% byte ref =[0 1 2 3 4 5]
% example = [34 0 0 164 31 0]
%
% example[5]*256*256*0.5+example[4]*256*0.5+example[3]*0.5+example[2]*0.5*1/256...
% example(4)*2^15+example(3)*2^7+example(2)*2^-1+example(1)*2^-9
description = cell(HDR.num_records,1); %24 bytes --> varies in size; some descriptions are longer, some are shorter
remainder = cell(HDR.num_records,1); % --> varies in size because description is not always the same length.
intro_size = 10; % 6 + 4
% intro_size = 34;
% remainder = zeros(HDR.num_records,bytes_per_record-intro_size,'uint8');
for r=1:HDR.num_records
start_sample(r) = fread(fid,1,'uint32'); % 4 bytes
evtID(r,:) = fread(fid,2,'uint8'); % [13,1]
padding = fread(fid,1,'uint32')'; %[0 0 0 0] % 4 bytes
curRecord = fread(fid,(bytes_per_record-intro_size)/bytes_per_uint16,'uint16')';
descriptionStop = find(curRecord==0,1,'first');
description{r} = char(curRecord(1:descriptionStop-1));
remainder{r}=curRecord(descriptionStop:end);
% description{r} = fread(fid,12,'uint16=>char')'; %24 bytes %need to read until I get to a 34 essentially%now at 64 or %32 bytes read
% remainder(r,:) = fread(fid,bytes_per_record-intro_size,'uint8')';
end
stop_sample = start_sample+5*embla_samplerate_out; %make it 5 seconds for viewing
elseif strcmpi(eventType,'custom')
%.evt is a 12 byte record
% 20 byte record
customType = stop_sample;
intro_size = 10;
if(false)
fseek(fid,32,'bof');
a=fread(fid,[bytes_per_record,inf],'uint8')';
fseek(fid,32,'bof');
a
end
% Earlier notes from Hyatt
%80 byte blocks
% 1 5 9 13 14 15 16 17 25 30 31 32 33 38 39 40 41 42 45 46 47 49 53 57 58 61 77 78 79 80
% [uint32][uint32][uint8*4][ uint8] [uint8][uint8 ][uint8] [uint8*8][uint8*5][uint8 ][uint8 ][uint8][uint8*5][uint8 ][uint8 ][uint8][uint8 ][uint8*3][uint8][uint8][uint8*2][uint32 ][uint32][uint8 ][uint8*3][uint32*4][uint8][uint8][uint8][uint8]
% [start ][stop ][7 2 3 0][1/32/33][ 0 ][0-248 ][0/63/64/191/192][ 255 ][ 0 ][0/64/128/192][0/86/87][0/64 ][0 ][0/64/128/192][0/84/85/86/87][ 0/64][ 44/46 ][ 0 ][0-255][0/1/2][ 0 ][? or 255][ 255 ][counter++/255][0/255 ][ 255 ][42 ][0/128][0-248][0/63/64/191/192]
% i
%intro_size = 8;