顺利更改日长

如何解决顺利更改日长

我想模拟一天长度随时间平滑变化（但保持正弦曲线）的样子。用于更改瞬时频率的“线性调频”的公式在https://en.wikipedia.org/wiki/Chirp中给出，但是在5天的24小时内编码，然后在另外5天的时间内过渡到12小时，则看起来不正确：>

period = list(  c(24,24,5),c(24,12,5) )
alpha = list(   c(0,c(0,5)  )
s_samples = 100
A=50
O=50
simulatedData = data.frame(t=numeric(),v=numeric()) #initialise the output
daySteps = c(0,cumsum(unlist(period)[seq(3,length(unlist(period)),by=3)])) #set up the period starts and ends to set over,starting at 0
##Cycle over each of the items in the list
for(set in seq(period) ){
  t_points = s_samples*period[[set]][3]
  t = seq(daySteps[set],daySteps[set+1],length.out=t_points) #make the time
  slope = (24/period[[set]][2]-24/period[[set]][1])/(max(t)-min(t)) # get the slope
  f0 = 24/period[[set]][1] - slope*(min(t)) # find the freq when t0
  c = (24/period[[set]][2]-f0)/(max(t)) #calculate the chirp see https://en.wikipedia.org/wiki/Chirp and https://dsp.stackexchange.com/questions/57904/chirp-after-t-seconds
  wt = ((c*(t^2))/2) + f0*(t) # calc the freq 
  a = alpha[[set]][1]
  v = A * cos(2*pi*wt - a) + O
  simulatedData = rbind(simulatedData,data.frame(t,v) )
}
plot(simulatedData,type="l",lwd=2)
t = seq(0,sum(unlist(period)[seq(3,by=3)]),by=1/24)
points(t,A*cos(2*pi*t)+O,col=3,lty=2)
points(t,A*cos(2*(24/12)*pi*t)+O,col=4,lty=2)

正如预期的那样，前24个是完美的，而后5天的最后一部分与周期为12h匹配，但是该时期的前一部分看起来相差180度。怎么了？

解决方法

我认为您正在使它变得比所需复杂得多。请记住，许多R函数已经被向量化。以下函数将在f0和f1之间的频率t0和t1之间产生线性chi，并带有一个可选的phi参数来指定您希望序列开始的周期：

chirp <- function(f0,f1,t0,t1,phi = 0,n_steps = 1000)
{
  C <- (f1 - f0)/(t1 - t0)
  x <- seq(t0,length.out = n_steps)
  y <- sin(2 * pi * (C / 2 * (x - t0)^2 + f0 * (x - t0)) + phi) # Ref Wikipedia
  data.frame(x,y)
}

当然，它还可以通过在两个相同的频率之间“ chi”来生成图表的静态前半部分，因此我们可以通过这样做获得图表上x，y点的数据框

df <- rbind(chirp(1,1,5),chirp(1,2,5,10))

这将导致：

plot(df$x,df$y,type = "l")

请注意，在5到10天之间有7.5个周期，因此，如果要平稳地继续频率2，则需要将phi参数设置为半周期（即pi）：

df <- rbind(df,chirp(2,10,15,phi = pi))

plot(df$x,type = "l")

请注意，如果线性调频脉冲发生在原始信号的偶数周期内，则线性调频脉冲信号和2 Hz信号的相位仅在n秒后匹配。对于奇数，相位将相差180度。这是线性chi的数学结果。要看到这一点，让我们使用函数在6秒内鸣叫，以便相位在10秒时匹配：

plot(df$x,type = "l")
lines(df2$x,df2$y,lty = 2,col = "green")
lines(df3$x,df3$y,col = "blue")
lines(df$x,df$y)