[R] need for help for solving operations in a vector

Makram Belhaj Fraj belhajfraj.makram at gmail.com
Tue Dec 29 12:33:25 CET 2015


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.makram at 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.pikal at 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.makram at 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.pikal at 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-bounces at r-project.org] On Behalf Of Makram
>> > Belhaj Fraj
>> > Sent: Wednesday, December 23, 2015 8:35 AM
>> > To: r-help at r-project.org; r-help-owner at 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]]



More information about the R-help mailing list