R Find the frequency and wavelength above a given value using the conditional value in data.table

MED is inserted below

Mre

date<-c('2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-01','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02','2001-01-02') time<-c('07:00:00 GMT','08:00:00 GMT','09:00:00 GMT','10:00:00 GMT','11:00:00 GMT','12:00:00 GMT','13:00:00 GMT','14:00:00 GMT','15:00:00 GMT','16:00:00 GMT','17:00:00 GMT', '18:00:00 GMT','19:00:00 GMT','20:00:00 GMT','21:00:00 GMT','22:00:00 GMT','23:00:00 GMT','00:00:00 GMT', '01:00:00 GMT','02:00:00 GMT','03:00:00 GMT','04:00:00 GMT','05:00:00 GMT','06:00:00 GMT','07:00:00 GMT','08:00:00 GMT','09:00:00 GMT','10:00:00 GMT','11:00:00 GMT','12:00:00 GMT','13:00:00 GMT','14:00:00 GMT','15:00:00 GMT','16:00:00 GMT','17:00:00 GMT','18:00:00 GMT','19:00:00 GMT','20:00:00 GMT','21:00:00 GMT') el<-c(0.257,0.687,1.861,3.288, 4.821,6.172,7.048,7.258,6.799,5.654,4.463,3.443,2.704,2.708,3.328,4.23,5.244,5.985,6.317,6.074,5.234,3.981,2.662,1.615,0.88,0.746,1.405,2.527,3.928,5.283,6.517,7.179,7.252,6.625,5.454,4.214,3.144,2.491,2.357) Time<-as.POSIXct(paste(date, time),tz="GMT") wave<-data.table(Time, el) ggplot(wave, aes(wave$Time, wave$el)) + geom_point() + labs(x="time", y="elevation") + geom_hline(aes(yintercept=4)) 

I have a time series of waves, and I want to have a function that can tell me the frequency and the average / average wavelength above a certain level. In my example, I chose 4.

I want to interpolate the time when the wave reaches 4 on the rising and falling edges and finds the time difference between two points for each wave.

I can do this with a for loop, but I think I could do it in data.table a lot faster. I have 1 mil + points for multiple locations and don't think the for loop will be efficient.

For the rising wave, I want to do something like:

 wave[,timeIs4:=ifelse(elev<3 & elev[+1]>4,TRUE,FALSE )] 

But instead of TRUE, enter my interpolation calculation. I do not know how to access the previous and outgoing values ​​in the data table, for example, in the for loop i + 1 or i-1.

Desired output

Rising foot I want to interpolate between points 4 and 5; 15 and 16; 29 and 30.

Falling leg I want to interpolate between points 11 and 12; 21 and 22; 36 and 37

Estimated result

 Rising Falling 10:28:00 17:27:00 21:45:00 3:59:00 11:03:00 18:12:00 

Then I can subtract Rising from Falling using difftime () to determine the time during which the water level was above a given level.

This will give me the frequency and duration of the water above a given level.

+5
source share
2 answers

Here a solution is possible using the devel version from GH . This will be needed for both the shift function (as indicated by @Jan) and the new dcast method, which accepts expressions. In addition, you do not have minutes in your MRE, so you do not know where you got them in the expected output.

In any case, for starters, we will create an index (we will call it Wave , so you know what wave it comes from #), which will tell us whether the wave grows or falls with shift . Then we will dcast according to the agreed values, removing unmatched with na.omit (you can change the order of the column names later if you like to use the setcolorder function)

 library(data.table) ## V 1.9.5+ dt[elev <= 4 & shift(elev, type = "lead") > 4, Wave := "Rising"] dt[elev > 4 & shift(elev, type = "lead") <= 4, Wave := "Falling"] dcast(na.omit(dt), cumsum(Wave == "Rising") ~ Wave, value.var = "time") # Wave Falling Rising # 1: 1 2001-01-01 17:00:00 2001-01-01 10:00:00 # 2: 2 2001-01-02 03:00:00 2001-01-01 21:00:00 # 3: 3 2001-01-02 18:00:00 2001-01-02 11:00:00 
+6
source

Here is another possible idea:

 elev = 4 #a helper function to calculate elapsed time ff = function(el1, el2, el, time1, time2) time1 + ((el - el1) / (el2 - el1)) * (time2 - time1) dif = diff(findInterval(wave$el, c(-Inf, elev, Inf))) ris = which(dif == 1) #risings fal = which(dif == -1) #fallings ff(wave$el[ris], wave$el[ris + 1], elev, wave$Time[ris], wave$Time[ris + 1]) #[1] "2001-01-01 10:27:52 GMT" "2001-01-01 21:44:42 GMT" "2001-01-02 11:03:11 GMT" ff(wave$el[fal], wave$el[fal + 1], elev, wave$Time[fal], wave$Time[fal + 1]) #[1] "2001-01-01 17:27:14 GMT" "2001-01-02 03:59:05 GMT" "2001-01-02 18:12:00 GMT" 
+3
source

All Articles