Hi Petr I runned the code you gave me as following, and I am adjusting for the threshold as suggested,
sss<-smooth.spline(qWC1, tint, nknots=length(tint)/2) peak <- which.max(predict(sss)$y) #qWC1=temp$theta mtime=temp$int baseline <- min(predict(sss)$y) vyska <- predict(sss)$y[peak] # 50% threshold HM <- (vyska+baseline)/2 plot(tint, qWC1,col=(tint>HM)+1, pch=19) #10% threshold HM<-vyska-(vyska*.1) plot(tint,qWC1, col=(tint>HM)+1, pch=19) cheers Makram On Tue, Dec 29, 2015 at 3:26 PM, Makram Belhaj Fraj < belhajfraj.mak...@gmail.com> wrote: > Hi Petr, > I apologize It was my first time using r-help so I didn't know how to > replay to email to all or not, > I am replaying to all for this email, > > Many thanks for the code, I am trying to use it, > please find attached the data file in csv, > > following are the first lines of code to read the data and calculate qWC1 > > sdata<-read.csv("almondc_10augv.csv",head=TRUE,sep=",") > > tint=sdata$scan #time intervall > > mtime=sdata$mtime #measurement time > > v1=sdata$vwc1 #value of moisture in percent > > qWC1=200*v1 #conversion in mm to get the variable I am working on > > > > On Tue, Dec 29, 2015 at 11:32 AM, PIKAL Petr <petr.pi...@precheza.cz> > wrote: > >> Hi >> >> >> >> You shall send your posts to the list, others can answer your questions >> and not only you can benefit from their answers too. >> >> >> >> As you did not post any data it is hard to say what are your issues. I >> believe that there are several values which are the same not only near the >> peak but also at the bottom part. If your data look like I remember and you >> want to keep all values near the peak value regardless they are slightly >> growing or falling, one approach can be to identify peak value, and select >> all values near the peak (the threshold is up to you). >> >> >> >> something like that >> >> >> >> sss<-smooth.spline(temp$theta, temp$int, nknots=length(temp$int)/2) >> >> peak <- which.max(predict(sss)$y) >> >> baseline <- min(predict(sss)$y) >> >> vyska <- predict(sss)$y[peak] >> >> # 50% threshold >> >> HM <- (vyska+baseline)/2 >> >> plot(temp$theta, temp$int, col=(temp$int>HM)+1, pch=19) >> >> #10% threshold >> >> HM<-vyska-(vyska*.1) >> >> plot(temp$theta, temp$int, col=(temp$int>HM)+1, pch=19) >> >> >> >> > dput(temp) >> >> structure(list(theta = c(28.995, 29.005, 29.015, 29.025, 29.035, >> >> 29.045, 29.055, 29.065, 29.075, 29.085, 29.095, 29.105, 29.115, >> >> 29.125, 29.135, 29.145, 29.155, 29.165, 29.175, 29.185, 29.195, >> >> 29.205, 29.215, 29.225, 29.235, 29.245, 29.255, 29.265, 29.275, >> >> 29.285, 29.295, 29.305, 29.315, 29.325, 29.335, 29.345, 29.355, >> >> 29.365, 29.375, 29.385, 29.395, 29.405, 29.415, 29.425, 29.435, >> >> 29.445, 29.455, 29.465, 29.475, 29.485, 29.495, 29.505, 29.515, >> >> 29.525, 29.535, 29.545, 29.555, 29.565, 29.575, 29.585, 29.595, >> >> 29.605, 29.615, 29.625, 29.635, 29.645, 29.655, 29.665, 29.675, >> >> 29.685, 29.695, 29.705, 29.715, 29.725, 29.735, 29.745, 29.755, >> >> 29.765, 29.775, 29.785, 29.795, 29.805, 29.815, 29.825, 29.835, >> >> 29.845, 29.855, 29.865, 29.875, 29.885, 29.895, 29.905, 29.915, >> >> 29.925, 29.935, 29.945, 29.955, 29.965, 29.975, 29.985, 29.995 >> >> ), int = c(329L, 330L, 318L, 287L, 315L, 344L, 333L, 324L, 334L, >> >> 366L, 339L, 374L, 375L, 335L, 415L, 371L, 413L, 382L, 408L, 406L, >> >> 407L, 440L, 475L, 465L, 516L, 510L, 490L, 550L, 663L, 647L, 628L, >> >> 721L, 789L, 814L, 890L, 923L, 1085L, 1102L, 1222L, 1356L, 1521L, >> >> 1729L, 1868L, 2120L, 2491L, 2656L, 3196L, 3599L, 4128L, 4536L, >> >> 5043L, 5310L, 5638L, 5792L, 5699L, 5374L, 4886L, 4473L, 4293L, >> >> 3757L, 3319L, 2934L, 2422L, 1998L, 1753L, 1397L, 1163L, 972L, >> >> 854L, 775L, 648L, 695L, 616L, 553L, 554L, 509L, 530L, 483L, 482L, >> >> 406L, 451L, 422L, 403L, 393L, 396L, 348L, 367L, 428L, 345L, 384L, >> >> 330L, 342L, 312L, 313L, 323L, 328L, 340L, 322L, 330L, 305L, 311L >> >> )), .Names = c("theta", "int"), row.names = 100:200, class = "data.frame") >> >> > >> >> >> >> Cheers >> >> Petr >> >> >> >> >> >> *From:* Makram Belhaj Fraj [mailto:belhajfraj.mak...@gmail.com] >> *Sent:* Tuesday, December 29, 2015 5:54 AM >> *To:* PIKAL Petr >> *Subject:* Re: [R] need for help for solving operations in a vector >> >> >> >> Hi Petr, >> >> I would like to thank you for your time >> >> i choose the following modification on the plotting so to keep only >> successive equals in red >> >> plot(qWC1, col=(c(0, diff(qWC1))==0 )+1) >> >> then, the red points I want to include in the irrigation period are the 3 >> successive red points in the plateau (equal values close to the max) >> >> These points are very important as they are corresponding to saturation, >> we continue to irrigate with the same flow and values remains constant for >> a certain time, so I must add them to irrigation >> >> >> >> up to now I didn't succeed in adding this plateau values to >> theirrigations using the diff(qWC1, lag=1) >> >> however, I wrote also some loops trying to catch all the irrigation and >> recharge separately and I still have some issues, >> >> following are the loops I used, with comments corresponding to the issues >> >> x<-qWC1 >> >> length(x) >> >> irrig<-rep(1,61) >> >> for (i in 2:61) { >> >> if (x[i-1]<x[i]){ >> >> irrig[i]<-x[i]-x[i-1] >> >> } >> >> } >> >> rech<-rep(1,61) >> >> for (i in 2:61) { >> >> if (x[i-1]>x[i]){ >> >> rech[i]<-x[i-1]-x[i] >> >> } >> >> } >> >> plot(x, type = "l", col = "black", ylim = c(min(0), max(92))) >> >> lines(irrig, type = "l", col = "cornflowerblue", ylim = c(min(0), >> max(15))) lines(rech, type = "l", col = "brown", ylim = c(min(0), max(15))) >> >> temp<-irrig<1.000001 #logical command to identify low values into a >> temporary vector temp >> >> temp2<-as.numeric(irrig>1.000001) #logical command to identify high >> values with 1 >> >> temp2 >> >> temp3<-as.numeric(rech>1.000001) >> >> temp3 >> >> irrig2<-irrig*temp2 #remove values inf 1 mm inirrig >> >> rech2<-rech*temp3 #same for rech >> >> plot(irrig2, type = "l", col = "cornflowerblue", ylim = c(min(0), >> max(15))) lines(rech2, type = "l", col = "brown", ylim = c(min(0), max(15))) >> >> Many thanks for your time and intellectual generosity >> >> Cheers >> >> Makram >> >> >> >> On Mon, Dec 28, 2015 at 12:15 PM, PIKAL Petr <petr.pi...@precheza.cz> >> wrote: >> >> Hi >> >> On top of answers you have got here is some plotting you need to answer >> yourself >> >> plot(qWC1, col=(c(0, diff(qWC1))>=0 )+1) >> >> Which from those red points you want to be included in irrigation period? >> All of them? Only part? Which part? >> >> Based on your figures you probably will not get 100% correct answer. >> >> Cheers >> Petr >> >> >> >> > -----Original Message----- >> > From: R-help [mailto:r-help-boun...@r-project.org] On Behalf Of Makram >> > Belhaj Fraj >> > Sent: Wednesday, December 23, 2015 8:35 AM >> > To: r-help@r-project.org; r-help-ow...@r-project.org >> > Subject: [R] need for help for solving operations in a vector >> > >> >> > Dear colleagues >> > i need your generous help to solve the following problem >> > >> > I have a soil moisture time series qWC1 (61 values) >> > > qWC1 >> > 75.33336 75.20617 75.20617 74.95275 74.95275 74.70059 74.70059 >> > 74.70059 >> > 74.57498 74.44968 74.32469 74.07563 85.57237 90.40123 90.73760 >> > 90.73760 90.73760 90.73760 90.90648 91.07582 91.24564 90.90648 86.82135 >> > 80.69793 >> > 79.30393 78.62058 78.21484 77.81226 77.67876 77.41279 77.28032 76.88495 >> > 76.75383 76.75383 76.49260 76.36249 76.23270 76.23270 76.10325 >> > 75.97412 >> > 75.84532 75.71685 75.71685 75.71685 75.71685 75.46087 75.46087 >> > 75.46087 >> > 75.33336 75.20617 75.20617 75.20617 75.20617 75.20617 75.20617 75.07930 >> > 75.07930 75.07930 74.95275 74.95275 74.95275 >> > >> > I want to measure consecutive increases corresponding to irrigation and >> > consecutive decreases corresponding to recharge I wrote the following >> > code and it does not calculate for each increment in i? >> > also note that I choose to not use diff command in time series because >> > I want also that "plateaux" corresponding to a minimum of 2 equal >> > consecutive values are accounted as positive differences=irrigations so >> > when x[i+1]==x[i] the difference y might be equal to the previous value >> > xi >> > >> > following the code i wrote >> > >> > x<-ts(qWC1,start=1, end=61, frequency=1) x[1] plot(x, type="h", col = >> > "green") >> > y<-rep(0,61) >> > for (i in 1:61) { >> > if (x[i+1] > x[i]){ >> > y[i]==x[i+1]-x[i] >> > } else if (x[i+1]==x[i]){ >> > y[i]=x[i+2]-x[i] >> > } else { >> > y[i]==x[i+1]-x[i] >> > } >> > >> > } >> > plot(y, type="h", col = "blueviolet") >> > >> > Many thank >> > Makram >> > >> >> >> >> ------------------------------ >> Tento e-mail a jakékoliv k němu připojené dokumenty jsou důvěrné a jsou >> určeny pouze jeho adresátům. >> Jestliže jste obdržel(a) tento e-mail omylem, informujte laskavě >> neprodleně jeho odesílatele. Obsah tohoto emailu i s přílohami a jeho kopie >> vymažte ze svého systému. >> Nejste-li zamýšleným adresátem tohoto emailu, nejste oprávněni tento >> email jakkoliv užívat, rozšiřovat, kopírovat či zveřejňovat. >> Odesílatel e-mailu neodpovídá za eventuální škodu způsobenou modifikacemi >> či zpožděním přenosu e-mailu. >> >> V případě, že je tento e-mail součástí obchodního jednání: >> - vyhrazuje si odesílatel právo ukončit kdykoliv jednání o uzavření >> smlouvy, a to z jakéhokoliv důvodu i bez uvedení důvodu. >> - a obsahuje-li nabídku, je adresát oprávněn nabídku bezodkladně >> přijmout; Odesílatel tohoto e-mailu (nabídky) vylučuje přijetí nabídky ze >> strany příjemce s dodatkem či odchylkou. >> - trvá odesílatel na tom, že příslušná smlouva je uzavřena teprve >> výslovným dosažením shody na všech jejích náležitostech. >> - odesílatel tohoto emailu informuje, že není oprávněn uzavírat za >> společnost žádné smlouvy s výjimkou případů, kdy k tomu byl písemně zmocněn >> nebo písemně pověřen a takové pověření nebo plná moc byly adresátovi tohoto >> emailu případně osobě, kterou adresát zastupuje, předloženy nebo jejich >> existence je adresátovi či osobě jím zastoupené známá. >> >> This e-mail and any documents attached to it may be confidential and are >> intended only for its intended recipients. >> If you received this e-mail by mistake, please immediately inform its >> sender. Delete the contents of this e-mail with all attachments and its >> copies from your system. >> If you are not the intended recipient of this e-mail, you are not >> authorized to use, disseminate, copy or disclose this e-mail in any manner. >> The sender of this e-mail shall not be liable for any possible damage >> caused by modifications of the e-mail or by delay with transfer of the >> email. >> >> In case that this e-mail forms part of business dealings: >> - the sender reserves the right to end negotiations about entering into a >> contract in any time, for any reason, and without stating any reasoning. >> - if the e-mail contains an offer, the recipient is entitled to >> immediately accept such offer; The sender of this e-mail (offer) excludes >> any acceptance of the offer on the part of the recipient containing any >> amendment or variation. >> - the sender insists on that the respective contract is concluded only >> upon an express mutual agreement on all its aspects. >> - the sender of this e-mail informs that he/she is not authorized to >> enter into any contracts on behalf of the company except for cases in which >> he/she is expressly authorized to do so in writing, and such authorization >> or power of attorney is submitted to the recipient or the person >> represented by the recipient, or the existence of such authorization is >> known to the recipient of the person represented by the recipient. >> > > [[alternative HTML version deleted]] ______________________________________________ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.